SATELLITE OBSERVATIONS REVEAL MORPHOLOGICAL CHANGES IN ASIAN DELTAS FROM 1986 TO 2018

: Asian deltas are densely populated and richly bio-diversified regions with significant social, economic, ecological and environmental importance across the globe. These low-lying and fertile floodplains support substantial food production. Because of sediment supply reductions, sea-level rising, subsidence, and frequent storm surges, many deltas are rapidly sinking and shrinking, which will threaten survival of millions of populations. It’s of great significance to monitor the delta morphological evolution and understand the relationship between their shoreline migrations and sediment supplies. However, quantifying the delta morphological changes is challenging due to its nature of high dynamics. In this study, we developed a methodology integrated by a time-series algorithm of water frequency and a spatial analysis to track the morphological changes of 11 major Asian deltas from 1986 to 2018.The results show that the deltas of Yellow, Yangtze, Pearl, Red, Mekong, Ganges-Brahmaputra, and Chao Phraya have experienced a net land gain at the total of 1829.3 km 2 , with China’s three major deltas accounting for 73 percent of all. The deltas of Mahanadi, Irrawaddy, Indus, and Godavari-Krishna have undergone a net land loss at the total of 295.1 km 2 , of which the Indus River Delta accounted for 65 percent. Most of the shorelines in the Yangtze River Delta were constantly advancing seaward at a rate of over 30 m/year, whereas, many ones in the rest of Asian deltas were significantly being eroded at a rate of over 20 m/year. These changes in delta morphology are closely related to sediment supplies and management policies. These results also suggest that rational wetland reclamation with auxiliary projects of ecological restoration can enhance the deltaic ability of sustainable development.


INTRODUCTION
River deltas are densely populated and richly bio-diversified regions with significant social, economic, ecological and environmental importance all over the world. In particular, Asia's mega-deltas support more than 400 million people and extensive biodiversity hotspots (Giosan et al., 2014). However, many deltas, such as the Mekong River Delta and the Indus River Delta, are rapidly sinking, shrinking, and facing numerous risks of increased rates of coastal erosion, frequent exposures to storm surges, heavy losses of vital habitats and associated ecosystems due to changing climates, sea-level rises and human activities (Besset et al., 2019;Minderhoud et al., 2019;Nienhuis et al., 2020;Syvitski et al., 2009). Among them, human activities have become a dominant factor in reshaping the deltaic morphology by regulating water discharges and sediment supplies (Warrick et al., 2019), destroying large-scale coastal wetlands (Tian et al., 2016), and developing estuarine navigable channels (Chen et al., 2021). Understanding the history of morphological evolution of these Asian deltas under such high-intensity of human activities presents significant implications for optimizing the managements and restorations of deltas in the rest of the world.
Healthy deltas are characterized by sediment depositions sufficient to keep up with the combined speeding of sea-level rise and surface subsidence (Anthony et al., 2014;Syvitski et al., 2009). That is, the fronts of deltas can continuously prograde, and the deltaic shorelines have been advancing seaward. Otherwise, the deltaic shorelines will be eroded and deltas will retreat. At present, satellite observations enable such processes to be captured with considerable accuracy and repeated monitoring.
Particularly, since the first launch of the earth observation satellite in 1972, Landsat archives have been considered as the most valuable data source for monitoring the dynamics of natural 1 resources and anthropogenic activities in both large spatial scales and long temporal scales (Masek et al., 2020). Based on the huge amount of data, deltaic morphodynamics have been investigated from different perspectives. Chen et al. (2020) provided the comprehensive remote sensing-based survey and review of the recent morphological evolution of the Irrawaddy River Delta under the anthropogenic impacts. Jarriel et al. (2020) identified the hotspots of morphological changes to the Ganges-Brahmaputra River Delta by calculating the channelized response variances on the basis of thirty-year Landsat imagery. Other morphological studies with remotely sensed imagery were also performed on the Mississippi River Delta (Ortiz et al., 2017), the Nile River Delta (Ghoneim et al., 2015), the Po River Delta (Ninfo et al., 2018), the Yangtze River Delta (Qiao et al., 2018), and the Mekong River Delta (Li et al., 2017) respectively.
Despite that remote sensing benefits the understanding of the deltaic evolution, existing studies above are limited by the methodological automation and applicability. Many approaches remain restricted to the manual delineation of deltaic shores from artificially selected cloud-free imageries. Therefore, they are time-consuming, subjective, and even unpracticable to detect long time-series of delta morphological changes. In addition, deltaic morphology is highly dynamic and easily affected by tides, waves, currents, and other anthropogenic disturbances. Detecting the delta morphological changes from a limited number of images may also yield inaccurate results. To overcome the influence of these mercurial factors on the description of deltaic tidal flats, some time-series-based approaches were well developed, such as the quantile synthesis method (Zhao et al., 2020) and the spectral index composite method . However, these approaches based on the time-series big data of earth observation are still rarely applied by the continental-scale studies of deltaic morphology.
In this paper, we aim to develop a rapid and robust method to quantify the morphological changes to the 11 Asian river deltas. Accordingly, the erosional or accretional areas and rates and spatial migrations of shorelines from 1986 to 2018 across the major deltas in Asia will be calculated respectively. The expected methodology can be adopted to study and assess the threatened global systems of deltas, from which the results shall provide profound insights into the developments, managements, and protections of other deltas in the rest of the world. The 11 Asian river deltas, including Yellow River Delta, Yangtze River Delta, Pearl River Delta, Red River Delta, Mekong River Delta, Mahanadi River Delta, Irrawaddy River Delta, Indus River Delta, Godavari-Krishna River Delta, Ganges-Brahmaputra River Delta, and Chao Phraya River Delta, are selected as the study cases (Figure 1(a)). These deltas, located in developing countries, are ones of the most densely populated regions in Asia and contain numerous mega-cities like Shanghai and Mumbai. They also support a great deal of agricultural developments, particularly for rice productions. It was estimated that Asia contributed about 90% of the global rice yields and the Mekong River Delta harvested over half of rice crops in Vietnam (Schneider and Asch, 2020). These Asian deltas differ in their dominating forces and divide into the river-dominated and tidedominated deltas. Recently, they are facing various degrees of threats from rapid decline in sediment supply, ground surface subsidence, and sea-level rise, and therefore exhibit different erosion-accretion patterns.

Data collection
In this study, we applied a total of 11,702 images of Landsat surface reflectance to analyze the morphological changes of 11 Asian deltas from 1986 to 2018 (Figure 1(b)). These images georeferenced with high accuracy were comprised of the Landsat 5 Thematic Mapper (TM), the Landsat 7 Enhanced Thematic Mapper-plus (ETM+), and the Landsat 8 Operational Land Imager (OLI). The Quality Assessment (QA) band of each Landsat image generated by the CFMask algorithm was adopted to identify and remove the ineffective observations (e.g., cloud, cloud shadows, etc.). Owing to the spatially and temporally consistent resolutions at the global scale of Landsat missions, a three-year temporal-scale water frequency composite imagery was created for each delta from the data covering seven periods of 1986-1988, 1991-1993, 1996-1998, 2001-2003, 2006-2008, 2011-2013, and 2016-2018. This is in order to minimize all the effects of cloud covers, tide fluctuations, and hydrological extremes (Chen et al., 2021). All the Landsat images were individually filtered and processed on the cloud-based Google Earth Engine (GEE) platform.

Quantification of delta morphological changes
To delineate the deltaic morphology from Landsat images and quantify their multi-temporal evolutions, a methodology framework integrated of time-series water index and spatial analysis was developed. This methodology consisted of four key steps ( Figure 2): (1) the generation of water body based on the water frequency algorithm; (2) the delineation of raster and vector of deltaic morphology; and (3) the calculation of delta morphological changes and shoreline migrations by using the Digital Shoreline Analysis System. Figure 2. Flowchart of the methodology adopted in this study.

Water frequency composite
Because of the shortwave infrared (SWIR) band being less sensitive to the sediment concentration and other optical active constituents within the water body than the near-infrared (NIR) band, the Modified Normalized Difference Water Index (MNDWI) would be, therefore, more stable and reliable than other water indexes in the estuarine areas (Huang et al., 2018). For each deltaic area, we calculated the time-series MNDWI from Landsat images during the related period. To distinguish water and land, threshold determination would be a critical issue.
To avoid the manual adjustments in this study, Otsu's method was employed to automatically calculate the threshold value of water/land segmentation from each MNDWI image. Water frequency images were composited to mitigate the impacts of tidal fluctuations and to extract the maximum areas of deltaic exposures. For time-series MNDWI of the related period in each delta, the water frequency at every pixel location was calculated as the quotient between the total frequency of the water pixel appearances at the position and the total frequency of observations at that position.

Deltaic morphology based on raster and vector
To delineate the deltaic morphology during low tide periods from the image sequences, 95% of threshold value was applied as the land-water segmentation, in other words, pixels with water frequency no less than 95% were considered as water bodies, otherwise, it would be regarded as potential land. The reason for not choosing a threshold value greater than 0 was to avoid certain errors, e.g., the extraction error of water body caused by cloud removal and cloud shadow. Moreover, the connected component analysis was applied in the resultant binary images to remove the irrelevant inland water and small islands. Vector-based deltaic maps were automatically generated by the ArcScan tool from the binary images, and then the results of automatic vectorization were manually examined and the deltaic shorelines were eventually extracted. The automatic vectorization tool in ArcGIS could avoid subjective errors caused by the manual digitization of shorelines.

Calculation of delta morphological changes
To quantify the delta morphological changes, the Digital Shoreline Analysis System (DSAS), an addition of ArcGIS that was developed by the USGS, was employed to compute the multi-temporal shoreline migrations throughout the period of 1986-2018. The DSAS required the integration of shorelines and initial baselines into a geodatabase. Thus, we constructed the baselines from the extracted deltaic shorelines and tried to keep two of them as close as possible. To calculate the shoreline movements, the DSAS cast the transect lines perpendicular to the baselines at defined intervals (200 m in this study). All of the transects were intersected with the multi-temporal shorelines while the associated intersections were applied as measurement points for calculating the rates of shoreline movements. The End Point Rate, as one output in the DSAS, was also calculated through dividing the movements by the time elapsed between the oldest and the newest shorelines, which therefore can be recognized as the shoreline migration rate. The net area changes of the 11 Asian deltas from 1986 to 2018, as well as the change during six periods, were calculated respectively ( Figure 3)  i.e., 6 km 2 /year. It should also be noted that the Pearl River Delta and the Chao Phraya River Delta both continually propagated seaward; however, the rates of their net land gain had already slowed down. Figure. 4 depicts the shoreline migrations of Asian deltas from 1986 to 2018. All of them showed more or less spatial and temporal variabilities in the migrations. Most of the shorelines across the Yangtze River Delta advanced seaward at a rate of over 30 m/year, and only the south bank of the Chongming Dongtan was slightly eroded and receded. Numerous shorelines along these deltas of the Yellow, Red, Mekong, Mahanadi, Indus, Irrawaddy, Godavari-Krishna, Ganges-Brahmaputra, and Chao Phraya rivers have experienced massive erosion at the rate of over 20 m/year. Major erosion was observed around the tidal inlets of these deltas. It should be note worthy that since the Yellow River changed its course, severe erosion has occurred near the old estuary and rapid propagation seaward has emerged around the new one. In this study, a three-year composite image for each period was generated to increase the accuracy of mapping deltaic morphology at the low tide. Despite the use of full time-series Landsat, the maximum portion of the extracted deltaic morphology may inevitably be smaller than that exposed during the actual low tide. Previous studies suggested that tidal correction would improve the accuracy of mapping deltaic morphology (Chen and Chang, 2009;Vos et al., 2019), but the accuracy of tidal data itself can also bring great uncertainty to the correction results. In addition, the 95% water frequency threshold was selected to automatically delineate deltaic morphology at periods of low tide. The choice of water/land segmentation threshold determines the spatial extents of deltas. Influenced by the dynamics of the tides, figure 5(a) shows the transition zone between absolute land and sea. When selecting a small threshold (e.g., 50%), part of the inter-tidal zone will be lost in the extraction results ( Figure 5(c)). While selecting a large threshold has no significant effect on the extraction of effective range of delta, but it will lead to a lot of noise errors ( Figure 5(e)), which will cause substantial inconveniences for the automatic processing. Thus, 95% water frequency was determined as the water/land segmentation threshold ( Figure 5(d)). To assess the accuracy of the 95% threshold for obtaining the delta shoreline, we compared it with the manually digitalized shoreline during low tide and found that the average shoreline displacement was no more than 30 m (i.e., the size of a Landsat pixel) ( Figure 5(b)).

Reasons causing spatiotemporal variations of delta morphology
The major reason behind deltaic coast erosion is the sediment trappings caused by the dams built at the upstream of rivers. It was estimated that there were more than 60,000 large dams worldwide, with over 3,700 more being currently planned or under construction. The sediment trappings from these dams presented a volume equivalent to an area of about 7,300 km 2 , assuming a 10m thick bed of sediments (Syvitski and Kettner, 2011). Sediment supply decreased more rapidly in the large Asian deltas, previous studies indicated that there had been ~90% decline in Yellow (Yu et al., 2013), ~70% decline in Yangtze (Yang et al., 2015), ~61% decline in Red (Ve et al., 2021), ~ 52% decline in Pearl (Zhang et al., 2011), ~ 74% decline in Mekong (Van Binh et al., 2020), ~67% decline in Mahanadi, ~94% decline in Indus, ~30% reduction in Ganga-Brahmaputra, ~87% decline in Krishna, ~74% decline in Godavari respectively (Gupta et al., 2012), over the past several decades. Currently, numerical simulation suggests that river sediment loads will continue to decline under the effects of anthropogenic activities and climate changes (Dunn et al., 2019). Such a situation may result in more deltaic coasts being under the threat of erosion, and trigger profoundly morphological, ecological, and biological responses in the deltas. The Yangtze River Delta, however, continues to advance seaward at the rate of over 30 m/year despite a plunge in its sediment load. Figure 6 shows the continuous accretions in the east shoal of Chongming Island, Yangtze estuary. In the past 35 years, the shoal has been reclaimed for four times. During each reclamation, dikes were built to separate the tidal environments and many natural wetlands were preserved outside the dikes. In this way, the preserved wetland vegetation still can dissipate waves and slow down flows, trap sediments and protect the marsh platform from erosion. Additionally, many projects of ecological restoration, such as the tidal flat filling with the soils dredged from deep-water channels, and the vegetation planting, were implemented in the tidal wetlands outside the dikes. These measures enable deltas to continually expand despite the sediment declines. Moreover, the reclaimed area can satisfy the rising land demands under rapid economic development. Different from the utilization pattern of the Yangtze River Delta, the inter-tidal wetlands of some Asian deltas were directly developed into the semi-open aquaculture ponds (exchanging waters with the oceans through the tidal cycles). For example, extensive mangroves were artificially removed and the tidal flats along the front of the Mekong River Delta were reclaimed into the aquaculture ponds, resulting in massive erosion along the south coasts of Mekong. By contrast, the pattern of "Reclamation-Restoration-Siltation" in Yangtze River Delta is therefore a more sustainable way to better balance the development and protection of river deltas.

CONCLUSION
In this study, more than three decades of spatiotemporal changes of morphologies in 11 Asian major deltas were individually studied by a water frequency algorithm and a spatial analysis method. The spatial extents and area variations of delta morphology, and the delineated shorelines at different periods were all quantitatively analyzed. The results demonstrated that these delta changes varied widely between a total net land gain of 1829.3 km 2 in 7 deltas and a total net land loss of 295.1 km 2 in 4 deltas. Numerous shorelines in Asian deltas were being significantly eroded at the rate of over 20 m/year due to the rapid declines in sediment supplies. Further drastic decreases of sediment load in the Asian deltas may trigger more profound impacts on their coasts and human communities in the future. Besides that, the differences in deltaic utilization and management policies can be directly reflected by the patterns of erosion and accretion to the deltas. The particular pattern of "Reclamation-Restoration-Siltation" in Yangtze will provide an excellent example of balancing the development and protection of river deltas.