Tracking the dynamics of tidal wetlands with time-series satellite images in the Yangtze River Estuary, China

ABSTRACT

Tidal wetlands provide a variety of ecosystem services to coastal communities but suffer severe losses due to anthropogenic activities in the Yangtze River Estuary (YRE). However, the detailed dynamics of tidal wetlands have not been well studied with sufficient spatiotemporal resolution. Here, we proposed a rapid classification method that integrates the COntinuous monitoring of Land Disturbance (COLD) algorithm and Median Composite (MC) based on the dense Landsat time series to track the dynamic processes of tidal wetlands in the YRE from 1990 to 2020. The results showed that the COLD-MC demonstrated remarkable effectiveness in detecting the change of tidal wetlands and excellent overall accuracy and kappa coefficient ranging from 90% to 96% and 0.89–0.95, respectively. The overall accuracy of change detection was 97% with an absolute error of 0.4 years. We found that the total area of tidal wetlands experienced a net loss of 59.75 km² in the YRE, but the gain and loss of the study period were 1556.07 and 1615.82 km², respectively. Land reclamation, sediment reduction, and Spartina alterniflora invasion pose significant threats to tidal wetlands. Sustainable management could be implemented through the establishment of nature reserves and ecological sediment enhancement engineering projects.

KEYWORDS:

• Tidal wetlands
• landsat time-series
• COLD-MC
• dynamic equilibrium
• land reclamation
• sediment starvation

1. Introduction

Tidal wetlands, thriving at the interface of terrestrial and marine ecosystems, primarily consist of mangroves, salt marshes, and mudflats (McLeod et al. 2011). They rank among the most productive ecosystems on the planet (Barbier et al. 2011; Davidson et al. 2019). Additionally, tidal wetlands play a crucial role in delivering a diverse range of irreplaceable ecosystem services, such as habitat provision for fauna and flora, biodiversity maintenance, coastal protection, erosion control, water purification, carbon sequestration, as well as raw materials and food provision, contributing to a staggering economic value of 20.4 trillion USD per year (Barbier et al. 2011; Castagno et al. 2021; Li et al. 2018; McLeod et al. 2011). Despite the vital economic and ecological importance of tidal wetlands, they are facing unprecedented pressures from climate change and anthropogenic activities (Murray et al. 2022). Over the period from 1999 to 2019, tidal wetlands experienced a net loss of 4,000 km² globally (Murray et al. 2022), resulting in a substantial decline in the ecological functions that they provide (Li et al. 2018). Anthropogenic activities are considered the primary factor contributing to the loss of intertidal wetlands, particularly in developing nations with rapid economic development. For example, just from 1975 to 2015, a substantial decline of 2,953 km² was observed in the tidal wetlands along coastal regions of China due to extensive land reclamation (Song et al. 2020). Except for direct occupation by anthropogenic influences, tidal wetlands are subject to the impacts of various other factors such as dam construction, invasion of alien species, and sea level rise (Liu et al. 2018; Ma et al. 2014; Wang, Li, and Xiang 2022; Yang et al. 2006). Under the combined influence of various factors (Newton et al. 2020), tidal wetlands have become highly vulnerable ecosystems, casting a veil of uncertainty over their future fate. Therefore, comprehending the causes, processes, and future of tidal wetland changes is of great significance for the conservation, management, and restoration of tidal wetlands.

Due to the capability of satellites to make long-distance and large-scale earth observations with short-term revisits, mapping tidal wetlands using remote sensing images has become the primary method to study the dynamic changes of tidal wetlands (Huang et al. 2022; Sun, Liu, et al. 2021; Zhu 2017). With more than half a century of archived earth observation data, dense Landsat time-series data (DLTS) provides the longest temporal records of earth observation data with considerable temporal and spatial resolution, showing great potential in improving our knowledge of the dynamics of tidal wetlands under global change from vast coverage and long-term scale (Wu, Zhi, et al. 2022; Wulder et al. 2022). Moreover, the emergence of cloud computing platforms, such as Google Earth Engine (GEE), has created an environment for the invocation and computation of DLTS (Amani et al. 2020). Tidal wetland mapping categories that incorporate statistical methods to reduce tidal inundation, vegetation phonology, and cloud cover are widely developed. For example, Jia et al. calculated the maximum values of different spectral indices to extract the range of tidal wetlands and separated tidal flats using the Otsu algorithm (Jia et al. 2021). Wang and Wu et al. counted the annual frequency of marine water and tidal vegetation based on all available Landsat images and generated the annual maps of tidal wetlands (Wang et al. 2020; Wu, Zhi, et al. 2022). These studies account for spectral information over a certain period, typically 1–3 years, to eliminate the influence of short-term tidal dynamics, resulting in accurate maps of tidal wetlands. However, converting the temporal information of pixels into statistical data ignores the temporal continuity and correlation of pixels, which may lead to misclassification of tidal wetland categories under the influence of tidal inundation and vegetation phenology (Di Vittorio and Georgakakos 2018; Sagar et al. 2018). Furthermore, the thematic mapping of tidal wetlands still faces the challenge of not fully harnessing the potential of remote sensing data which record the entire developmental trajectory of tidal wetlands (Yang et al. 2022).

In recent years, several change detection algorithms have been developed to detect land use and land Cover (LULC) changes to make the most comprehensive use of all available remote sensing data. Such algorithms, including Landsat-based detection of Trends in Disturbance and Recovery (Landtrendr) (Kennedy, Yang, and Cohen 2010), Breaks For Additive Seasonal and Trend (BFAST) (Verbesselt et al. 2010), Continuous Change Detection and Classification (Zhu and Woodcock 2014), and COntinuous monitoring of Land Disturbance (COLD) (Zhu et al. 2020), detect changes by creating time-series models to track the development trajectories of each pixel (Zhu 2017). Based on the single band time series data, the Landtrendr and BFAST can detect external disturbances on a pixel scale (Rhif et al. 2022; Wu, Zhi, et al. 2022). Originally designed to monitor disturbances in terrestrial forests, Landtrendr and BFAST have found expanded utility in the monitoring of changes within intertidal wetlands (Chen et al. 2022; Wu, Zhi, et al. 2022). In comparison to Landtrendr and BFAST, CCDC and COLD, which calibrated based on various disturbance events including wetlands, combine multiple spectral bands to capture more LULC changes with higher accuracy (Zhu et al. 2020; Zhu and Woodcock 2014). The divergence lies in the fact that the COLD algorithm has updated the model initialization and change detection rules within the framework of the CCDC algorithm, leading to a notable enhancement in the precision of change detection (Zhu et al. 2020). In recent years, CCDC and COLD algorithm have been widely applied in the research of dynamic changes in tidal wetlands ecosystems, demonstrating substantial promise in unraveling the intricate dynamics of tidal wetland transformations (Awty-Carroll et al. 2019; 2021; Yang et al. 2022). However, condition change introduces the possibility of spurious changes being detected by the COLD algorithm, especially in active estuarine coastal regions, thereby potentially deviating the classification results from the ground truth (Yang et al. 2022). Moreover, the COLD algorithm has only been applied to a single category of mangrove or salt marsh vegetation, and it is less involved in mapping tidal wetlands of complicated classes. This issue underscores a critical challenge necessitating resolution within the scope of this study.

The Yangtze River Estuary (YRE) represents a pivotal ecological region with extensive tidal wetlands that serve as crucial supporters of multiple ecosystem services (Wu et al. 2003), and it stands as a prominent center of economic development in China (Sun, Zhu, and Cheng 2018). However, the rapid expansion of YRE's economy has attracted substantial human settlement and engendered profound anthropogenic stressors on tidal wetland habitats (Chen et al. 2018). Additionally, the Yangtze River Estuary is one of the areas that are most severely affected by the invasion of Spartina alterniflora in China (Liu et al. 2020). The integrated influences of human activities and climate changes have precipitated rapid transformations within the tidal wetlands of YRE (Tian et al. 2015). Consequently, delving into the dynamic fluctuations in the tidal wetlands of YRE serves as a representative marker, mirroring the broader panorama of tidal wetland changes not just within China but also on a global scale. However, due to the highly dynamic nature and high turbidity of the YRE, conventional mapping methods face challenges in generating a fine-scale tidal wetland spatial distribution dataset with high spatiotemporal resolution.

To address the aforementioned issues, we propose a rapid and robust classification method that can be applied to tidal wetlands with highly dynamic change and high turbidity based on statistical methods and change detection algorithm and generate a high spatiotemporal resolution spatial distribution dataset of tidal wetlands in the YRE covering the period from 1990 to 2020 using all available DLTS. From the results with higher spatiotemporal resolution derived from DLTS data, we aim to answer the following two questions:

1. How applicable is the classification method in highly dynamic estuarine coastal regions?

2. What are the distribution, dynamic changes, and drivers of the tidal wetlands in the YRE?

2. Data and methods

2.1. Study area

YRE is situated between Shanghai and Jiangsu Province, with geographical coordinates ranging from 30°33'38″ to 32°10'23″ N and 120°56'25″ to 122°16'15″ E. The Yangtze River carries abundant sediment into the estuary, with an annual average sediment discharge of 122 million tons in the last decade (Ministry of Water Resources of China 2022). YRE falls under the classification of a subtropical monsoon climate, characterized by an annual mean temperature of 17.3°C and precipitation of 1409.1 mm (Shanghai Statistical Bureau 2020). Ample sediment and hydrological conditions have contributed to the formation of diverse tidal wetland ecosystems, including Nanhui Dongtan, Jiuduansha, Chongming Dongtan, and Chongming Beitan. As habitats and wintering grounds for many migratory birds, Jiuduansha and Chongming Dongtang were listed as national nature reserves in 2005. The Yangtze River Delta is one of the fastest-growing economic regions in China. In recent years, as urbanization continues to accelerate, many large-scale reclamation activities and water conservancy projects have taken place along the coast, seriously squeezing the habitats of tidal wetlands. In this study, the land boundary was defined by the coastline in 1990, while the marine boundary was defined by the 15-km-wide buffer zone from the coastline in 2020. The area between the land boundary and the marine boundary is determined as the study area (Figure 1).

Figure 1. (a) The red rectangle represents the study area, the Yangtze River Estuary. (b) The spatial distribution of the number of good-quality observations used in the COLD algorithm. (c) The number of Landsat images from 1990 to 2020.

2.2. Data and processing

2.2.1. Pre-processing of the DTLS

In this study, all available Landsat collections from 1 January 1990 to 31 December 2020, including Landsat 5, Landsat 7, and Landsat 8 images with a total number of 3227 were used to explore the dynamics of tidal wetlands in the YRE. We obtained these data from the United States Geological Survey official website (https://glovis.usgs.gov/). All images were geo-registered and inter-calibrated high-quality data across the Landsat instruments, making it well-suited for time series analysis (Awty-Carroll et al. 2019). Blue, Green, Red, Near Infrared (NIR), Short-wave Infrared 1 (SWIR1) and Short-wave Infrared 2 (SWIR2) bands were used in this study.

2.2.2. Sample dataset

Firstly, proportionable samples were generated according to the area of different LULC types in the YRE based on the GEE (Peng et al. 2023; Stehman and Xing 2022). Secondly, a multi-source incorporating high-resolution images from Google Earth, Sentinel-2 images and relevant reference, was employed for visual interpretation of these samples. In the absence of high-resolution data from historical periods, monthly Landsat images were synthesized using GEE. The interpretation of samples was conducted by considering various band combinations, taking into account spectral features, phenological characteristics, and geographical information. A total of 14,309 sample were obtained, including 1921 Phragmites australis samples, 1339 Scirpus mariqueter samples, 1454 Spartina alterniflora samples, 1684 tidal flats samples, 2129 marine water samples, 781 impervious surface samples, 1085 aquaculture ponds samples, 1211 land water samples, 1295 farmland samples and 1410 other LULC samples. The spatial distribution of the sample points in 2020 is shown in Figure 2. Finally, we utilized the stratified sampling method to select validation and training samples (Jia et al. 2021).

Figure 2. Spatial Distribution of Sample Points in 2020. 1-Phragmites australis, 2-Scirpus mariqueter, 3-Spartina alterniflora, 4-tidal flats, 5-marine water, 6-impervious surface, 7- aquaculture ponds, 8-land water, 9-farmland, 10-other.

2.2.3. Auxiliary data

The shoreline data for the YRE from 1990 to 2020 were obtained by visually interpreting Landsat images taken during high tide periods, with a temporal resolution of 2 years (Wu, Gao, et al. 2022). The social-economic data, such as Gross Domestic Product (GDP) and the resident population during the same period were requested from the Shanghai Bureau of Statistics (https://tjj.sh.gov.cn/index.html). We also collected the sediment data for the Yangtze River during the years 2000–2020 were sourced from the Information Center of the Ministry of Water Resources in China (http://xxzx.mwr.gov.cn/xxgk/gbjb/zghlnsgb/).

2.3. Methods for extracting tidal wetlands

The detailed workflow is presented in Figure 3. First, time series models and synthesized annual spectral features were developed using the COLD algorithm and median composite method, respectively. Second, a random forest algorithm was employed to integrate the established time series models and synthesized annual spectral features as input features for the classifier, generating thematic maps of tidal wetlands from 1990 to 2020 in the YRE. Third, an analysis of the tidal wetland map was conducted based on the land use transition matrix and centroid migration model. Finally, the dynamic changes in the tidal wetlands at the YRE over the past three decades were delineated.

Figure 3. The workflow for this study.

2.3.1. Time-series models construction

In this study, the COLD algorithm, which was developed for continuous monitoring of land disturbance and long-term LULC classification for large areas (Awty-Carroll et al. 2021; Zhu et al. 2020), was used for constructing time-series models of tidal wetlands in the YRE over the period from 1990 to 2020. The COLD algorithm employs the Fourier model to create time series models for various bands on a pixel-by-pixel basis, and then compares new observations to predicted values to detect LULC changes. Each time series model (Equation (1)) consists of a trend term and three harmonic components, where the trend term could capture the inter-annual variation characteristics of the time series, and the harmonic components are used to estimate the intra-annual variation characteristics of the time series (Xie et al. 2022). (1) $\begin{array}{c}\hat{\rho }(i,x)=a_{0,i}+c_{1,i}x+\sum _{k=1}^3\left\{\left(a_{k,i}\cos {\displaystyle \frac{2\mathrm{\pi x}}{T}+b_{k,i}\sin {\displaystyle \frac{2\mathrm{\pi x}}{T}}}\right)\right\}\end{array}$(1) where $\hat{\rho }(i,x)$ is the predicted value of surface reflectance for the ith Landsat band at the xth Julian day; $a_{0,i}$ is the overall value of surface reflectance for the ith Landsat band (i.e. intercept); $c_{1,i}$ indicates the inter-annual change of the ith band (i.e. slope); $k$ represents the time frequency of harmonic components ($k$ = 1, 2, and 3); $a_{k,i}$ and $b_{k,i}$ denote the intra-annual change of the ith Landsat band; $T$ means the number of days per year, which is 365.25.

The operational procedure of the COLD algorithm is delineated into three principal steps. In the first step, stability testing is conducted during the model initialization phase because the changes occurring within this stage can impact the subsequent predictions and change detection of the model. Changes may occur at the beginning, middle, and end stages of model initialization, thereby causing alterations in the slope of the model. Model stability is evaluated based on Equation (2), which combines the absolute value of the slop of the time series model ($L_{\mathrm{slop}}$) and the larger absolute values between the difference in prediction value and clear Landsat observations at the start and end ($L_{\mathrm{start}}$ and $L_{\mathrm{end}}$). The optimal outcome is attained when the change probability is set to 0.99 (Zhu et al. 2020). If the model is unstable, COLD automatically incorporates next clear observation and removes the earliest clear observation until the model stabilizes. (2) $\begin{array}{c}|L_{\mathrm{slop}}|+max\{|L_{\mathrm{start}}|+|L_{\mathrm{end}}|\}<\mathrm{change}\mathrm{probability}\end{array}$(2) In the second step, following stability testing, COLD incorporates new clear Landsat observation to update the existing time series model and detect changes by comparing the differences between clear observation and predicted value. The distinction between clear observation and predicted value is determined in COLD based on Equation (3). COLD assumes that the sum of the squares of model deviations follows a chi-square distribution (${\chi }^2$), where the degrees of freedom are equal to the number of spectral bands used for change detection. Subsequently, a change probability of 0.99 is established based on the probability derived from the chi-square distribution (Equation (3)). Once six consecutive clear observations exceed the change probability threshold, COLD will flag it as a potential change. Given the approximate alignment of changes in direction, the COLD algorithm further confirms whether a real change has occurred by calculating the angle between each pair of neighboring change vectors for these six anomaly observations. The COLD algorithm specifies that a change event can be confirmed when the average angle of change vectors is less than 45° (Equation (4)). (3) $\begin{array}{c}\sum _{i=1}^k{\left({\displaystyle \frac{\rho (i,x)-\hat{\rho }(i,x)}{\mathrm{RMS}{E}_i}}\right)}^2\sim {\chi }^2>{\chi }_{0.99}^2\end{array}$(3) where $\rho (i,x)$ and $\hat{\rho }(i,x)$ is the actual predicted value of surface reflectance for the ith Landsat band at the xth Julian day; $k$ is the number of bands used for change detection. (4) $\begin{array}{c}{\displaystyle \frac{1}{k}\sum _{i=1}^{k-1}{\beta }_{i,i+1}<{45}^{\circ }}\end{array}$(4) where $i$ is the ith anomaly observation; $k$ means the number of consecutive anomaly observations. ${\beta }_{i,i+1}$ represents the neighboring change vectors between the ith and i + 1th anomaly observation.

Finally, once a change is detected, the time series model undergoes disruption, resulting in a breakpoint representing the occurrence of the change event. Afterward, a new time series model is established to search for the next breakpoint, and this process continues until clear Landsat observations are exhausted. All the spatial–temporal information of time series models on a pixel-by-pixel basis will be recorded in a specific data format and processed into raster files (Table 1).

Table 1. Output of the COLD algorithm.

Parameters	Description
t_start	The time when the model starts
t_end	The time when the model ends
t_break	The time when the change event occurred
coefs	Coefficients for each time series model for each spectral band, including $a_{0,i}$, $c_{1,i}$, $a_{k,i}$ and $b_{k,i}$ (Equation (1)).
rmse	Root mean square error for the model fit for each spectral band
change_prob	Change probability (0–1). Change probability equals to the number of consecutive anomaly observations divided by N (N: number of consecutive anomaly observations to confirm a change)
pos	Position of each time series model in the image (x,y)
num_obs	Number of available clear observations used for constructing model

An instance of change detection utilizing the COLD algorithm for a pixel (31.7553°N,121.59°E) which undergone LULC change once over the past three decades in Figure 4. The COLD algorithm has constructed two segments of the time series model spanning the periods 1997–2018 and 2018–2020 for this pixel across multiple spectral bands, capturing the LULC change that occurred in 2018 (Figure 4(a)). By examining historical high-resolution images from Google Earth Pro, it is revealed that the LULC type of the pixel was identified as marine water from 1990 to 2017. However, in 2018, it underwent reclamation and transformed into land (Figure 4(b)). Table 1 describes the output variables of the COLD algorithm.

Figure 4. Illustration of the COLD algorithm. A pixel (31.7553°N,121.59°E) has undergone LULC change once from marine water to land over the past three decades. (a) LULC changes identified through the COLD algorithm. (b) LULC changes observed on high-resolution images from Google Earth Pro.

2.3.2. Synthetic annual images

We noticed the presence of errors within the constructed time series model (Figure 5). These deviations, arising from dynamic tidal inundation and non-LULC type alterations, introduced condition changes (Yang et al. 2022). The accumulation of these errors will lead to a decline in the accuracy of thematic maps. Therefore, we attempted to reflect the genuine ground conditions of each year by synthesizing annual images and integrating them into the classifier to mitigate such errors. Currently, the mainstream approach for the annual composite of remote sensing images in tidal wetlands is the maximum value composite (MVC) method, known for effectively delineating the boundaries of tidal wetlands (Jia et al. 2021). However, in estuaries with high turbidity like the YRE, the use of MVC amplifies the spectral information of turbidity (Figure 6(a)), potentially leading to an over-identification of tidal flats in classification (Jia et al. 2021). Conversely, the median synthesis method weakens the influence of turbidity by considering the annual median value, reducing the risk of misidentification (Figure 6(b)). Moreover, we found that MVC-derived images frequently exhibit noise and blank regions and images synthesized using the median composite (MC) method exhibit higher quality with minimal noise and blank areas (Figure 6). The MC sets each pixel in the composite based on which pixel in the time series Landsat image collection has the median value of spectral bands. In this study, the synthesized annual images consist of six bands, including SWIR1, SWIR2, NIR, RED, GREEN, and BLUE, which almost maintain consistency with the bands used for constructing the time-series model and encompass most of the spectral information from the ground (Zhu et al. 2020).

Figure 5. Validation of the detected change year. A comparison between the detected change years and the interpreted change years for sampling points. The numbers inside the squares represent the number of sample points. MAE: mean absolute error.

Figure 6. Annual images synthesized through different methods. (a) Maximum value composite (MVC). (b) Median composite (MC).

2.3.3. Classification

We utilized the Random Forest (RF) algorithm to classify each time series model into 10 LULC types from 1990 to 2020 using Google Earth Engine (GEE) which is a cloud platform that enables rapid classification and efficient processing of large volumes of data (Amani et al. 2020). RF is considered one of the most effective machine-learning classifiers for LULC classification due to its ability to adapt to diverse datasets and account for correlations and interactions among features (Chen and Ishwaran 2012). The number of decision trees in the random forest is set to 100, with mtry defaulting to the square root of the variable count. The LULC types were categorized into tidal wetlands (Phragmites australis, Scirpus mariqueter, Spartina alterniflora, and tidal flats) and none-tidal wetlands (impervious surface, aquaculture ponds, land water, farmland, marine water, and others (unused land and terrestrial vegetation)) (Table 2). The input features for the classifier consist of all coefficients and root mean square error (RMSE) of the time series model (Table 1) and annual images synthesized in section 2.3.2. The LULC of the tidal zone was classified from 1990 to 2020 with an interval of two years. Since the model was initialized in the initial year, the accuracy of classification may be affected. Due to the minimal changes observed between the satellite images of 1991 and 1990, we used the image of 1991 to represent the classification of 1990.

Table 2. Tidal wetland categories in the YRE.

Category I	Category II	Description
Tidal wetlands	Phragmites australis	salt marsh vegetation on high tidal flats
	Spartina alterniflora	salt marsh vegetation on middle tidal flats
	Scirpus mariqueter	salt marsh vegetation on low tidal flats
	Tidal flats	Bare mudflats without vegetation
None-tidal wetlands	Impervious surface	Buildings, roads, etc.
	Aquaculture ponds	Man-made ponds for fisheries production
	Farmland	Land used for agricultural production
	Land water	Inland rivers, lakes, and reservoirs
	Marine water	Seawater near the estuary
	Others	Terrestrial vegetation and unused land

2.4. Accuracy assessment

The accuracy assessment strategy for change maps and LULC maps is based on the stratified random sampling approach (Olofsson et al. 2014). For the change map, a total number of 600 random sampling points were generated, consisting of 387 stable points and 213 changed points. Additionally, 4209 random sampling points were generated for LULC maps. Then, all sampling points were interpreted with very high-resolution (VHR) imageries accessible in Google Earth Pro. Based on stratified reference samples, confusion matrixes are calculated including the user’s accuracy (UA), producer’s accuracy (PA), kappa coefficient (Kappa), and overall accuracy (OA) for change map and LULC maps (Olofsson et al. 2014). We also calculated the F1 score for the change map and the average of PA and UA (APA and AUA) for LULC maps based on confusion matrixes. Furthermore, a comparative analysis was conducted between the detected change years and the interpreted change years, aiming to provide a comprehensive evaluation of the accuracy of the change results.

2.5. Dynamic analysis method for tidal wetlands

2.5.1. Land use transition matrix

The land use transition matrix reveals both the direction and quantity of transitions between different LULC types, offering valuable insights into the spatiotemporal dynamics inherent in land use patterns (Eq (5)). In this study, the land use transition matrix is employed to investigate the spatial distribution of gains and losses in tidal wetlands of the YRE, quantifying the transitions into and out of tidal wetlands during different periods. (5) $\begin{array}{c}{S}_{\mathrm{ij}}=\left[\begin{array}{cccc}{S}_{11}& S_{12}& \cdots & S_{1n}\\ S_{21}& S_{22}& \cdots & S_{2n}\\ ⋮& ⋮& ⋮& ⋮\\ S_{n1}& S_{n2}& \cdots & S_{\mathrm{nn}}\end{array}\right]\end{array}$(5) where $S$ represents the area of land use; n represents the number of land use types; $i$ and $j$ represent the land-use types before and after the transition; $S_{\mathrm{ij}}$ represents the area of land use type $i$ transitioning to land use type $j$ during the study period.

2.5.2. Centroid migration model

The Centroid migration model can evaluate the spatiotemporal dynamics of tidal wetlands by computing centroids across different years and quantifying inter-centroid distances (Eq (6)&(7)). The centroid serves as a representation of the spatial distribution of tidal wetlands, while the migration distance of centroids serves as a quantitative measure of the extent of change within tidal wetlands. Based on the mean center tool in ArcGIS 10.7, we compute the centroid migration of tidal wetlands in the YRE from 1990 to 2020. (6) $\begin{array}{c}{\mathrm{X}}_{\mathrm{t}}={\displaystyle \frac{{\sum }_{i=1}^n{M}_{\mathrm{ti}}\times X_{\mathrm{ti}}}{{\sum }_{i=1}^n{M}_{\mathrm{ti}}}}\end{array}$(6) (7) $\begin{array}{c}{\mathrm{Y}}_{\mathrm{t}}={\displaystyle \frac{{\sum }_{i=1}^n{M}_{\mathrm{ti}}\times Y_{\mathrm{ti}}}{{\sum }_{i=1}^n{M}_{\mathrm{ti}}}}\end{array}$(7) where $X_t$ and $Y_t$ represent the longitude and latitude coordinates of the centroid of tidal wetland area distribution in year $t$, respectively; n denotes the number of landscape patches of tidal wetland. $M_{\mathrm{ti}}$ represents the area of the ith patch in year $t$; $X_{\mathrm{ti}}$ and $Y_{\mathrm{ti}}$ represent the longitude and latitude coordinates of the ith patch in year $t$, respectively.

3. Results

3.1. Accuracy assessment

3.1.1. Performance of the COLD algorithm in tidal wetlands

All pixels were categorized into stable and changed pixels, and the results of change detection were validated based on a confusion matrix (Table 3). The Overall accuracy of the change map is 97%, with a kappa coefficient of 0.94. The PA and UA for stable and changed pixels are all above 94%. The F1 scores reached 0.96 and 0.98. Results also showed that the detected change years exhibited a high degree of consistency with the interpreted change years, suggesting the robust capability of the COLD algorithm in accurately capturing the occurrence time of majority change events (Figure 5). The majority of deviations from the actual occurrence of change events fall within the range of 1–3 years, and the mean average error in validating the change detection results is 0.4 years. The application of the COLD algorithm for change detection demonstrated remarkable outcomes.

Table 3. Confusion matrix of change detection results.

	Stable pixel	Changed pixel	Total number	UA
Stable pixel	209	12	221	97.21%
Change pixel	6	373	379	96.88%
Total number	215	385	600
PA	94.57%	98.42%
F1 score	0.96	0.98
	OA = 97.00%		Kappa = 0.94

PA: producer’s accuracy, UA: user’s accuracy, OA: Overall accuracy, Kappa: kappa coefficient.

3.1.2. Accuracy assessment for LULC maps

Figure 7 shows the accuracy assessment of LULC maps from 1990 to 2020. The OA of LULC maps ranged from 90% to 96% and the Kappa ranged from 0.89 to 0.95. The APA and AUA of classification maps also achieved satisfactory results with high values of around 90%. With the exception of specific years, the consistency in precision evaluation metrics across different years remains relatively stable on the LULC maps. The results demonstrated that the LULC map is highly consistent with real ground conditions, which can be used for subsequent analyses.

Figure 7. The accuracy assessment of classification maps. OA: Overall accuracy, APA: Average producer accuracy, AUA: Average user accuracy.

3.2. Status and historical change of tidal wetlands in YRE

Tidal wetlands exhibit extensive distribution across YRE, encompassing vast expanses of wetland habitats in prominent regions including Chongming Beitan (Figure 8(a)), Chongming Dongtan (Figure 8(b)), Jiuduansha (Figure 8(c)), and Nanhui Dongtan (Figure 8(d)). The total extent of tidal wetlands is 667.5 km² in 2020, including 95.06 km² of Phragmites australis, 80.9 km² of Scirpus mariqueter, 151.1 km² of Spartina alterniflora, and 340.44 km² of tidal flats, respectively. Both Scirpus mariqueter and tidal flats are distributed throughout the YRE. Phragmites australis populations predominantly inhabit the southern expanse of Chongming Dongtan and Shansha of Jiuduansha, while Spartina alterniflora is mainly distributed in Chongming Beitan, and Zhongsha and Xiasha of Jiuduansha, as well as Nanhui Dongtan.

Figure 8. Type and spatial distribution of tidal wetland habitats in YRE in 2020. (a) Chongming Beitan. (b) Chongming Dongtan. (c) Jiuduansha. (d) Nanhui Dongtan. (e) Proportions of different LULC categories in YRE. The colors in e correspond to the legend.

The total area of tidal wetlands in YRE decreased slightly from 727.25 km² to 667.5 km² from 1990 to 2020. However, it exhibited distinct patterns over the past three decades, which can be divided into three stages: (1) gradual expansion from 1990 to 1996, (2) rapid decline from 1996 to 2010, and (3) a period of accelerated growth from 2010 to 2020 (Figure 9). During the initial stage, the tidal wetlands showed a trend of increasing by 10.64 km²/a. In the subsequent stage, the wetlands experienced a remarkable decline with a slope twice as steep as the previous stage (Figure 9), resulting in a loss of 383.06 km² of tidal wetlands with an average rate of 23.55 km²/a. Notably, this phase was characterized by the emergence of the invasive species of tidal wetland vegetation, Spartina alterniflora, encompassing an area of 54.17 km². The final stage witnessed a rebound of tidal wetlands, with tidal wetlands expanding at a rate of 18.67 km²/a. The overall vegetation area of tidal wetlands increased from 292.92 km² to 327.09 km² over the period from 1990 to 2020. However, the invasive species, Spartina alterniflora, gradually became the dominant wetland vegetation with an area of 151.1 km² in 2020, almost approaching the total extent of Phragmites australis and Scirpus mariqueter with an area of 175.96 km². Although there have been modest recoveries in the areas of Phragmites australis, Scirpus mariqueter, and tidal flats, they have yet to reach their original extents. Particularly noteworthy is the diminished presence of Scirpus mariqueter, which now occupies less than half of its extent in 1990.

Figure 9. Changes in tidal wetlands from 1990 to 2020.

The tidal wetlands also exhibited seaward movement. Notably, the largest movement occurred in Nanhui Dongtan, and the seaward movement of the centroid was 8.54 km from 1990 to 2020 (Figure 10(d)), while Jiuduansha is the smallest with a distance of 1.66 km (Figure 10(c)). Moreover, Chongming Beitan (Figure 10(a)) and Chongming Dongtan (Figure 10(b)) also showed seaward movement with the distance of 2.57 and 3.43 km, respectively. Consequently, over the last three decades, there has been a remarkable reduction in marine water, reaching an extensive area of 853.16 km².

Figure 10. Spatial displacement of area-weighted centroids for main tidal wetlands in YRE from 1990 to 2020. (a) Chongming Beitan. (b) Chongming Dongtan. (c) Jiuduansha. (d) Nanhui Dongtan.

3.3. Gains and losses of tidal wetlands in YRE

The gains and losses of coastal wetlands fluctuated in the YRE (Figure 11). Over the period spanning 1990–2020, an extensive area of 1556.07 km² was gained, while a substantial extent of 1615.82 km² was concurrently lost, ultimately resulting in a net loss of 59.75 km². The total ratio of wetland gains to losses is 0.96, highlighting a quasi-equilibrium that existed in the dynamics of tidal wetlands in the YRE over the study period. The tidal wetlands in YRE suffered the most severe losses between 1998 and 2000 and experienced the greatest gains from 2010 to 2012 with a total area of 158.81 and 170.01 km², respectively.

Figure 11. (a) Gains and losses of tidal wetlands from 1990 to 2020. (b) Gain/loss ratio of tidal wetlands from 1990 to 2020.

The gains and losses of tidal wetlands exhibit a spatial variation (Figures 12 and 13). The gains of tidal wetlands due to the natural deposition were characterized by the gradual expansive belt towards the sea, such as Chongming Island and Jiuduan Island (Figure 12). However, the gains of coastal wetlands were little in the position where along the main channel of Yangtze Estuary, such as the Southen Chongming Dongtan Wetland near to the North Channel, and Southern Jiuduansha along the South Passage (Figure 12(b and c)). The significant gains of coastal wetlands from 2016 to 2020 in Nanhui Dongtan Wetland (Figure 12(d)) resulted from the sedimentation promotion project. Almost 104.35 km² of tidal wetlands were formed rapidly.

Figure 12. Gain of tidal wetlands in YRE during 1990–2020. (a) Chongming Beitan. (b) Chongming Dongtan. (c) Jiuduansha. (d) Nanhui Dongtan.

Figure 13. Loss of tidal wetlands in YRE during 1990–2020. (a) Chongming Beitan. (b) Chongming Dongtan. (c) Jiuduansha. (d) Nanhui Dongtan.

The losses of tidal wetlands are affected by extensive land reclamation and natural erosion. A net area of 587.42 km² of coastal wetlands in the YRE was reclaimed from 1990 to 2020, which most showed regular shapes in space (Figure 13). Moreover, a substantial decline in tidal wetland area was observed in spatial, with a notable emphasis on Chongming Dongtan (Figure 13(b)) and Nanhui Dongtan (Figure 13(d)). The prevailing trend of area loss was predominantly directed from land to sea, primarily occurring before 2010. Additionally, a continuous loss of area from sea to land due to erosion was also identified in the southeastern region of Xiasha of Jiuduan Sha (Figure 13(c)).

The inflow and outflow of tidal wetlands have undergone various changes in the past 30 years (Figure 14). During the periods spanning from 1990 to 2000 and from 2000 to 2010, a declining trend in tidal wetlands was observed. Interestingly, a contrasting trend of wetland expansion emerged during the period from 2010 to 2020. During the period from 1990 to 2010, a considerable number of tidal wetlands were reclaimed, with respective areas of 274.66 and 205.41 km². A noticeable decline in reclamation activities occurred after 2010, leading to a reduced extent of tidal wetland reclamation, encompassing an area of 107.35 km². Among all the reclaimed tidal wetlands, ‘Others’ (unused land and terrestrial vegetation) represents the main type of wetland conversion. Over the period spanning 1990–2010, this category witnessed large conversion extents of 127.99 and 105.49 km², which declined to 47.12 km² during 2010–2020. Notably, the conversion to ‘Others’ accounted for approximately 50% of the total reclaimed area in each period. The remaining reclaimed areas consist of impervious surfaces, aquaculture ponds, land water, and farmland. Marine water is also one of the major types of tidal wetlands conversion, with areas of 85.33 and 92.94 km² being converted during the periods of 1990–2010. After 2010, the conversion of marine water decreased significantly to 24.17 km². Moreover, marine water is the primary type contributing to the influx of tidal wetlands. Over the three periods, it replenished tidal wetlands with areas of 209.55, 139.36, and 374.73 km², respectively. The contribution of other LULC types to the influx of tidal wetlands can be deemed inconsequential.

Figure 14. LULC transition in the tidal zone from 1990 to 2020. Tidal wetlands: Phragmites australis, Scirpus mariqueter, Spartina alterniflora, and tidal flats. MW: marine water; IS: impervious surface; AP: aquaculture ponds; LW: land water; FL: farmland; OT: others (unused land and terrestrial vegetation). The features marked with black dashed rectangular boxes are tidal wetlands. The thickness of the lines represents the area transferred in and out.

4. Discussion

4.1. Adaptability of the method

4.1.1. Reliability of tidal wetland mapping

In this investigation, we presented a method named COLD-MC for mapping the tidal wetlands in YRE based on all available Landsat images with highly accurate. The successful implementation of this investigation can be attributed to open data sources, robust algorithms, and the powerful computational capabilities of the GEE platform. Firstly, the availability of open-source Landsat time-series data has provided a wealth of clear observational data, offering the potential for continuous monitoring of highly dynamic tidal wetlands. Secondly, the COLD algorithm maximizes the use of all available clear pixels by continuously tracking the developmental changes of each pixel, provides temporal information about tidal wetlands, and effectively mitigates the short-term impacts caused by dynamic tidal inundation. Concurrently, owing to the elevated sediment concentration in estuaries, marine water reflectance experiences an increase in specific bands at certain times within a specific year. The MC algorithm weakens the influence of high turbidity in estuaries and avoids the occurrence of extreme values by selecting median values compared to MVC (Table 4 & Figure 15(a)), generating high-quality images representing the annual ground conditions and significantly reducing classification errors attributed to the time-series models constructed by the COLD algorithm (Table 4 & Figure 15(b)). For example, under the same conditions, MC provides a more comprehensive extraction of the marginal tidal flats (Figure 15(a)), achieving a more thorough separation between turbid marine water and tidal flats (Figure 15). The integration of continuous change detection algorithms and statistical methods effectively mitigated the uncertainty arising from environmental factors (Figure 15 COLD_MC). Lastly, the GEE platform offers an environment for extensive remote sensing image access, processing, and computation (Amani et al. 2020). The approach demonstrated remarkable success within the challenging zones of the YRE which are characterized by high turbidity and dynamic environmental shifts, paving the way for potential applications in broader regions in the future.

Figure 15. Classification results of different strategies in Jiuduansha. The distinctions in classification are denoted within the confines of the black rectangular frame. The reference Landsat 8 image was acquired on 2020/05/12. R:G:B = 5:4:3.

Note: The results here have not undergone any post-processing.

Table 4. Different strategies for classification accuracy.

	COLD	COLD-MVC	COLD-MC
Overall accuracy	88.00%	88.79%	92.80%
Kappa coefficient	0.84	0.85	0.91

4.1.2. Comparison with other datasets

We compared the tidal wetland vegetation and tidal flats in our tidal wetland map (TWM) with two other publicly available datasets, respectively (Figures 16–18). One dataset was obtained from Wang et al., who calculated the annual frequency of marine water and tidal wetland based on all available Landsat images and applied a decision tree to generate the annual map of tidal wetlands in YRE in 2020 at 30-m spatial resolution (FUDAN/OU) (Wang et al. 2023). Another dataset was obtained from Jia et al., who extracted a 10-meter resolution map of tidal flats of China in 2020 using the maximum composite method and OUST algorithm (CTF) based on Sentinel-2 images between 1 January 2019 and 30 June 2020 (Jia et al. 2021).

Figure 16. Comparison of tidal wetland areas in the YRE between our tidal wetland map (TWM) in 2020 and the FUDAN/OU map in 2020 (Wang et al. 2023) and the CTF map in 2020 (Jia et al. 2021). (a) Comparison of vegetation area in tidal wetlands. (b) Comparison of tidal flats area in tidal wetlands

Figure 17. Subset views of tidal wetland vegetation from our tidal wetland map (TWM) and the FUDAN/OU map in 2020 (Wang et al. 2023) in Chongming Beitan (A), Jiuduansha (B) and Nanhui Dongtan (C). The background Landsat 8 image was acquired on 2020/08/16. R:G:B = 5:4:3.

Figure 18. Subset views of tidal flats from our tidal wetland map (TWM), the FUDAN/OU map in 2020 (Wang et al. 2023), and the CTF map (Jia et al. 2021) in Chongming Beitan (A), Jiuduansha (B) and Nanhui Dongtan (C). The background Landsat 8 image was acquired on 2020/05/12. R:G:B = 5:4:3.

In comparison to the FUDAN/OU map, both datasets exhibited a remarkably high consistency in tidal wetland vegetation (Figure 16(a)), as FUDAN/OU reported the area to be 350 km², and the area we extracted was 327.06 km². Additionally, their spatial distributions exhibited a striking congruence (Figure 17). The slight differences in tidal wetland vegetation areas between the TWM and FUDAN/OU datasets may be attributed to mapping methods. The tidal wetland vegetation in FUDAN/OU was separated through a decision tree. In the YRE, Scirpus mariqueter and tidal flats, which are located on low tidal zones, tend to intermingle. Consequently, utilizing a decision tree often leads to misclassifying these mixed pixels (Wang et al. 2020). However, the FUDAN/OU dataset (200.39 km²) detected fewer tidal flats than our study (340.44 km²) (Figure 16(b)). The primary factor accounting for this disparity is attributed to mapping methods. Due to the influence of clouds and fog in coastal areas, clear Landsat pixels are already scarce (Zhi, Wu, and Su 2022). Moreover, the impact of tidal inundation may result in a higher frequency of pixels being classified as water rather than tidal flats (Figure 18 B_FUDAN/OU, C_FUDAN/OU). Therefore, under the same data source conditions, eliminating short-term tidal inundation effects necessitates considering the long-term developmental trajectory of pixels (Figure 18).

The TWM dataset had a much smaller tidal flats area than the CTF dataset, with 340.44 and 437.73 km², respectively (Figure 16(b)). The differences between the TWM and CTF datasets can be attributed to the data sources and mapping methods. Despite the higher spatial and temporal resolution of Sentinel-2, the practice of compressing short-term data into a group of maximum values data may result in the loss of essential temporal information in tidal wetlands (Yang et al. 2022), such as the phenological information of vegetation. A considerable portion of tidal wetland vegetation in the CTF dataset had been categorized as tidal flats, leading to an overestimation of the areas (Figure 18).

4.1.3. Uncertainties in tidal wetland mapping

The uncertainties of tidal wetland maps were primarily caused by image data, mapping algorithms and inappropriate training data. Firstly, despite utilizing the QA band to eliminate low-quality pixels, it is impossible to attain entirely high-quality pixels, which may introduce uncertainty into the resultant maps (Wang et al. 2020). Additionally, the quantity of clear observations for different pixels ranges from 114 to 705 (Figure 1), with an inconsistent spatial distribution that may impact the accuracy of change detection algorithms. Multiple data sources will be considered to enhance clear observations in future research, such as dense Sentinel-2 images. Secondly, the median composite method may lead to an underestimation of tidal flats, which could contribute to the smaller area of tidal flats in the TWM compared to the CTF (Figure 16(b)). Thirdly, given the analogous features of Phragmites australis and Spartina alterniflora in remote sensing images, visual interpretation of sample points is solely reliant on phenological characteristics. In regions characterized by mixed vegetation, the potential for misclassification may contribute to a decline in overall classification accuracy (Liu et al. 2022). Moreover, pixel-based classification is susceptible to pronounced salt-and-pepper artifacts. A label similarity probability filter may be used for post-classification to improve accuracy (Sun, Yang, et al. 2021). Finally, despite utilizing all available time-series remote sensing images to accurately extract intertidal flats, the duration of local tidal flats covered by marine water is significantly higher than their exposed duration, resulting in incomplete extraction of some local tidal flats (Jia et al. 2021; Peng et al. 2024; Zhang, Xu, et al. 2022).

4.2. Challenges faced by tidal wetlands under rapid dynamic changes

Intensive land reclamation has been conducted due to accelerated urbanization in the YRE, which resulted in substantial losses of tidal wetlands (Cheng et al. 2023; Sun, Zhu, and Cheng 2018; Tian et al. 2016). A significant area of 587.42 km² of tidal wetlands has been reclaimed for urban expansion (Figure 19(I)), agricultural cultivation (Figure 19(II)), aquaculture (Figure 19(III)), and the construction of reservoirs (Figure 19(IV)), which accounts for the decrease of intertidal wetlands in the second phase of Figure 9. In addition, reclamation-induced loss of intertidal wetlands is nearly irreversible, with only 5.39 km² of reclaimed land reverting to tidal wetlands (Figure 14) (Liu et al. 2020; Tian et al. 2015; Wu et al. 2018). Consequently, following extensive land reclamation, the majority of tidal wetlands undergo significant displacement from their original positions (Figure 10).

Figure 19. Typical land reclamation of the tidal zone in YRE. (I) Pudong Airpot; (II) Qianshao farm; (III) aquaculture pond; (IV) Qingcaosha Reservoir.

Declining sediment supply caused by large-scale dam construction along the Yangtze River is also one of the reasons for changes in intertidal wetlands, particularly following the operation of the Three Gorges Dam (TGD) (Yang et al. 2021). During the period from the initial impoundment of the TGD in 2003 to its peak water level in 2010 (Wang et al. 2013), the sediment supply experienced a pronounced decline (Figure 20), slowing down the expansion rate of tidal wetlands and even resulting in loss of subaqueous delta in the YRE (Yang et al. 2015). After the TGD reached its maximum water level, the sediment supply to the YRE peaked in 2012, possibly contributing to an increase in the area of tidal flats that year (Figure 9). During the period of reduced sediment, the vegetation within tidal wetlands has exhibited a sustained and relatively stable rate of accretion, particularly with an elevation rise in the vertical dimension (Shi et al. 2022; Yang et al. 2008; Zhang, Zhang, et al. 2022). In contrast, tidal flats suffer from a rapid decline in sediment rates, with localized areas transitioning from accretion to erosion due to the sediment rate decreases to a level lower than the rate of sea-level rise (Figure 13(c)) (Crosby et al. 2016; Shi et al. 2022). This transition is particularly pronounced in the tidal wetland zones between −5-meter and 0-meter isobath, where a discernible reduction in the area has been observed (Shi et al. 2022). Therefore, the surface of tidal wetland was progressively becoming steeper and narrower, impeding outward expansion through sediment deposition.

Figure 20. Sediment flux at Datong station and the area of tidal wetlands in the YRE.

The invasion of Spartina alterniflora is also recognized as one of the key factors driving dynamic changes in tidal wetlands (Liu et al. 2018). In 1997, the Shanghai municipal government introduced a limited area of Spartina alterniflora to Jiuduansha through ecological engineering, and it was not until 2000 that the presence of Spartina alterniflora was observed from Landsat images. Spartina alterniflora expands from 3.86 km² in 2000–151.1 km² in 2020 at a rate of 7.36 km²/a, and constitutes 46.2% of the tidal wetland vegetation in the YRE, making it the dominant vegetation in the tidal zone (Figures 8 and 9). Therefore, Spartina alterniflora plays a crucial role in maintaining the dynamic equilibrium of tidal wetland areas. However, the rapid expansion of Spartina alterniflora comes at the expense of displacing native vegetation habitats, especially Scirpus mariqueter (Shi et al. 2022). The Scirpus mariqueter community provides critical habitats and food for diverse waterbird species (Ma et al. 2017). As the habitat of Scirpus mariqueter communities diminishes, there has been a discernible reduction in the abundance and diversity of waterbird species (Ge et al. 2009; Ma et al. 2009; Zou et al. 2016). In addition, the prominent silt-promoting capacity of Spartina alterniflora leads to a rapid increase in elevation, while tidal creeks within the community exhibit increased density and reduced connectivity (Hou et al. 2013; Zong et al. 2023). Furthermore, the biogenic barriers formed by Spartina alterniflora communities restrict the inward extension of tidal creeks, hindering the entry of tidal water into Spartina alterniflora salt marshes or higher tidal flats (Hou et al. 2013; Shen, Zhang, and Wang 2003; Zong et al. 2023). Consequently, a cross-shore circulation pattern emerges on the outer side, triggering sediment erosion at the marsh edge (Shen, Zhang, and Wang 2003).

4.3. Protection and restoration of tidal wetlands

Despite the sustained loss of tidal wetlands throughout the past three decades, the net losses of tidal wetlands in the YRE amount to a mere 59.75 km², accounting for less than one-tenth of the extent in 1990. This is attributed to the efforts of the local government in the protection and restoration of tidal wetlands, such as establishing nature reserves and implementing ecological sediment promotion engineering projects (ESPEP) (Yuan et al. 2022). Jiuduansha Wetland National Nature Reserve and ESPEP in Nanhui Dongtan wetland serve as typical representatives. In 1997, the government introduced Spartina alterniflora and Phragmites australis in Jiuduansha. Subsequently, the establishment of the Jiuduansha Wetland National Nature Reserve has provided favorable conditions for the natural expansion of wetlands in Jiuduansha while effectively minimizing anthropogenic disturbances. During the past three decades, Jiuduansha wetland has experienced notable natural expansion, witnessing a substantial increase in its extent from 38.96 km² in 1990–155.33 km² in 2020 at a rate of 3.88 km²/a (Figure 21). Contrasting the protracted and gradual natural progression of wetland development, ESPEP has expedited the process of wetland formation. The construction of groins and breakwaters in ESPEP has effectively weakened hydrodynamic forces (Wu, Yang, et al. 2022), facilitating sediment accumulation and accelerating the formation of tidal wetlands within the project area (Yuan et al. 2022). Nanhui Dongtan wetland, characterized as one of the most severely reclaimed areas in YRE, has formed rapidly after the conduction of the ESPEP (2013–2018), with a growth rate of 13.38 km²/a, almost elevenfold compared to the area observed in 2012 (Figure 22). Notably, Jiuduansha and Nanhui Dongtan wetlands alone contribute to 40% of the total tidal wetland area in the YRE. However, a considerable abundance of Spartina alterniflora was observed in both Jiuduansha Wetland National Nature Reserve and Nanhui Dongtan wetland. Therefore, the proposition of Spartina alterniflora removal within the YRE and the implementation of the native wetland vegetation restoration project are under consideration.

Figure 21. The dynamic changes of Jiuduansha Wetland National Nature Reserve from 1990 to 2020. (a) Landsat image of Jiuduansha wetland. (b) Bar chart illustrating the changes in the area of Jiuduansha wetland. Subfigures (c)∼(f) represent the spatial changes in the area of Jiuduansha wetland in 1990, 2000, 2010, and 2020, respectively.

Figure 22. The dynamic changes of Nanhui Dongtan wetland from 2012 to 2020 (a) The breakwaters, groins, and shoreline of Nanhui Dongtan wetland. (b) Bar chart illustrating the changes in the area of Nanhui Dongtan wetland. Subfigures(c)∼(f) represent the spatial changes in the area of Jiuduansha wetland in 2012, 2014, 2018, and 2020, respectively.

5. Conclusion

In this study, we developed a rapid classification method for tidal wetlands, integrating COLD algorithms with MC, and utilized all available Landsat remote sensing images to track the spatiotemporal changes and explore driving factors for the tidal wetlands of the YRE over the past 30 years. The results showed that our method achieved satisfactory classification results in the tidal wetlands with an accuracy range of 90% to 96%. 667.5 km² of tidal wetlands were present in the tidal zone of YRE in 2020, a decrease of 59.75 km² compared to 1990. Anthropogenic activities were the primary drivers, which resulted in the migration of the tidal wetlands habitat towards the sea. The refinement of relevant legal and policy frameworks has led to a reduction in the intensity of land reclamation activities. The diminished sediment supply and the invasion of Spartina alterniflora have also exerted substantial influences on the genesis and ecological structure of tidal wetlands, respectively. Concurrently, the construction of natural reserves and the implementation of ESPEP have resulted in the formation of new tidal wetland habitats, mitigating a significant portion of the previously incurred losses. The discoveries in this study provide crucial insights for the government and stakeholders, enabling them to tackle the multifaceted challenges associated with the management, restoration, conservation, and sustainable development of the tidal wetlands in the YRE.

Acknowledgment

This research was funded by the National Natural Science Foundation of China (Grant No. 42301540), the Open Research Fund of State Key Laboratory of Estuarine and Coastal Research (Grant number SKLEC-KF202307), the Natural Science Foundation of Fujian Province (Grant No. 2022J05024), the Education Department of Fujian Province (Grant No. JAT210027), and the Natural Science Foundation for Distinguished Young Scholars of Fujian Province, China (Grant No. 2021J06014).

Data availability statement

The data presented in this study are available on request from the corresponding author.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This research was funded by the National Natural Science Foundation of China (grant number 42301540), the Open Research Fund of State Key Laboratory of Estuarine and Coastal Research (grant number SKLEC-KF202307), the Natural Science Foundation of Fujian Province (grant number 2022J05024), the Education Department of Fujian Province (grant number JAT210027), and the Natural Science Foundation for Distinguished Young Scholars of Fujian Province, China (grant number 2021J06014).