Tidal flats (TFs) are dynamic, coastal areas located within the transition zone, composed of mangroves, salt marsh, and mudflats, where they are regularly exposed and submerged by tidal fluctuations (Wang et al., 2020). Coastal TFs include mudflats, rocks, and sandy beaches (Dyer et al., 2000), serving as transitional zones that exhibit marine and terrestrial environments (Murray et al., 2019; Wei et al., 2015; Yao, 2013). TFs perform essential ecological functions, including pollutant absorption, shoreline stabilization, storm buffering, and the conservation of biodiversity (Buczko et al., 2022; Zhang et al., 2022).
However, TFs, which are vital to both ecological and economic functions, are increasingly threatened by human activities such as land reclamation and natural disturbances, such as global warming-induced sea-level rise (Kolker et al., 2011; Rodríguez et al., 2017). Anthropogenic activities, such as dam construction, industrial development, and aquaculture, have increased the vulnerability of TFs (Chen et al., 2020). Therefore, accurate and efficient monitoring of the spatial and temporal distribution of TF ecosystems is essential for assessing their ecological health and ensuring the protection of coastal zones. Accurately mapping TFs presents unique challenges due to their dynamic nature, as they are only fully exposed for brief periods during low tide and their characteristics rapidly change with tidal fluctuations (Zhao et al., 2022).
Early TF monitoring usually depended on field surveys using traditional equipment to map spatial extents (Zhao et al., 2020), but these methods were often unable to capture dynamic changes over time due to their labor-intensive nature and large coverage requirements. Advancements in remote sensing (RS) have significantly enhanced image acquisition, data storage, time-series observation capabilities, and computational efficiency (Campbell & Wang, 2019), which has improved TF mapping at the national and international scale (Murray et al., 2019; Wu et al., 2023). High-resolution photography and cutting-edge technologies, such as topobathymetric LiDAR and synthetic aperture radar (SAR) (Bell et al., 2016), are frequently used in local-scale TF. However, the application of these techniques to national-scale mapping remains limited due to performance constraints and insufficient data availability. TF mapping methodologies employing optical RS data can be broadly classified into two primary computational approaches: (1) supervised machine learning techniques and (2) rule-based decision tree algorithms leveraging domain-specific spectral and phenological thresholds.
TF categorization has been done using training sample-based machine learning techniques, including random forest (Zhang et al., 2019) and support vector machine algorithms (Zhou et al., 2013). However, these techniques require a large number of training samples. Additionally, they typically result in misclassifications of inland features, such as aquaculture ponds and mudflats (Murray et al., 2019). Knowledge-based decision tree methods have been used to classify TFs (Zhao et al., 2023), relying on empirical knowledge and expert-defined thresholds. Large-scale mapping approaches (e.g., Catalão & Nico, 2017; Australian mapping by Sagar et al., 2017 and Murray et al., 2019) have been developed. Among these, Zhao et al. (2020, 2023) introduced a quantile composite method, which utilizes the 95th and 5th percentiles from time-series imagery to capture the highest and lowest tidal conditions at each pixel. However, although effective for defining minimum tidal coverage, the quantile method underestimates the actual maximum spatial extent of TFs at high tide. This approach was used in a more recent study (Jia et al., 2021; Zhang et al., 2022), except that they used maximum rather than quantiles to create composite images. While this method performs well overall, it still faces challenges of spectral confusion.
The primary challenge in RS monitoring of TFs is the mismatch between the instantaneous nature of satellite imagery and the dynamic characteristics of TFs, leading to inaccuracies in resulting maps. However, mapping of TFs effectiveness is still limited by the temporal resolution of available satellite imagery and spectral similarities between TFs and neighboring wetland types (Gong et al., 2019). Nonetheless, the use of freely available satellite imagery and cloud-based platforms, such as Google Earth Engine (GEE), has greatly improved monitoring capabilities by enabling large-scale and facilitating historical land cover assessments across extensive spatial domains. In light of these limitations, a long-term TF monitoring, this study aims: (1) to develop a novel pixel-based composite approach utilizing time-series Landsat (1984–2024) and Sentinel-2 (2015–2024) imagery processed through GEE to quantify spatiotemporal dynamics of TF ecosystems in the coastal regions of China and Pakistan; (2) combines quantile composite, maximum spectral index composite (MSIC), and the Otsu Algorithm into a unified workflow, enables high-resolution in small, environmentally diverse coastal zones without relying on elevation data; (3) applies the proposed method to investigate the influence of adjacent land use and urbanization on TFs over the past 40 years; and (4) evaluates the accuracy of TF extraction by comparing the results with existing TF datasets.
Both sites are economically important and close urban centers with active coastal development, making them ideal for studying human impacts on TF dynamics. Although not directly compared, the two regions differ in tidal regimes and anthropogenic pressure, making them a representative basis for evaluating the technique across coastal settings. In addition to exposing the long-term impacts of urbanization on the coast, this research produces an improved TF map for a variety of coastal habitats and offers recommendations for strategies for coastal development and biodiversity preservation.
This study focuses on analyzing the dynamics of TFs in two geographically and environmentally distinct coastal regions in China and Pakistan. Both economically important sites are near major urban centers with active coastal development, making them relevant for understanding human activity and TF dynamics. The first study area is located in Pakistan along the eastern portion of the Karachi coast bordering the northern Arabian Sea in Fig. 1A. It extends to the port and east of Karachi Harbor (Chinna Creek), governed by the Karachi Port Trust (KPT) in Karachi, Pakistan. Several industrial districts and the main commercial district are close to its location, which is between the towns of Keamari and Saddar. It spreads between latitude 24°46′N and longitude 66°59′E, falls under a semi-arid climatic zone (Ahsanullah et al., 2021), and receives an average of 217 mm of rainfall annually. The tidal regime in this area is classified as mesotidal, exhibiting average tidal ranges of 1.8 m to 2.5 m, shaped by semidiurnal and mixed tidal patterns (https://hydrography.paknavy.gov.pk/tide-table/). The country’s major maritime hub, Karachi Port, handles 25.45 million tons of dry cargo and 11.74 million tons of liquid cargo yearly, accounting for 60% of commerce.
Geographical locations of study area. (A) KPT, Karachi, Pakistan. (B) Yangma Island, Yantai China. KPT, Karachi Port Trust.
The secondary study area comprises Yangma Island (37°23′0″N–121°36′0″E), a 13.52 km2 oval-shaped island located off Yantai′s coast in northeastern Shandong Peninsula, China (Fig. 1B). The shoreline measures roughly 19.5 km and exemplifies typical temperate TF ecosystems characterized by mixed sand and stony bottoms. It experiences a humid subtropical climate with distinct seasons, moderate rainfall, and high humidity (Mou & Liu, 2025) with approximately 524.9 mm annual average precipitation. The tidal range in this region is comparatively modest, typically categorized as microtidal, with average ranges of approximately 1.2 m. Over the past two decades, it has experienced considerable shoreline transformation attributed to the expansion of aquaculture, tourism development, and coastal engineering activities (Yang et al., 2022).
Distribution of sample points. True color composite image of Landsat for both study areas.
Landsat, a series of satellites deployed by the National Aeronautics and Space Administration (NASA), has been operational since 1972 (Crawford et al., 2023). These satellites provide moderate temporal (16 days) and spatial (30 m) resolution, accumulating extensive and long-term surface RS imagery that is commonly applied to Earth observation research. This research utilized Landsat images from 1984 to 2024, including thematic mapper (TM) and operational land imager (OLI) data acquired from the GEE platform. A total of 705 images for the KPT region in Pakistan and 502 images for Yangma Island in China were selected based on <20% cloud cover. All data were L2T level and masked for clouds, shadows, and snow using the quality assessment (QA) band to ensure temporal consistency. A data gap exists from 2003 to 2007 for KPT, as shown in Fig. S1 in Supplementary Materials, limiting the availability of high-quality imagery during this period. All valid images within a 3-year window were aggregated to construct representative composites, ensuring sufficient coverage of both high and low-tide states.
The Sentinel-2 mission comprises twin polar-orbiting satellites (Sentinel-2A launched June 2015; Sentinel-2B launched March 2017) operating in a sun-synchronous orbit. Equipped with the multispectral instrument (MSI), the constellation acquires data across 13 spectral bands (440–2200 nm) at three spatial resolutions: 10 m (visible and NIR: bands 2, 3, 4, and 8), 20 m (red edge/SWIR: bands 5, 6, 7, 8A, 11, and 12), and 60 m (atmospheric correction bands). It observes worldwide land surface changes with a large swath and a 2–5-days revisit cycle (Wang & Atkinson, 2018). Around 2500 Sentinel-2 images were selected with minimal cloud cover (0%–20%) from 2017 to 2024. The QA bitmask band was used to mask out bad-quality observations caused by opaque and cirrus clouds in each image to ensure a good number of quality images.
The sample points were manually gathered from Landsat and high-resolution Imagery from Google Earth, covering categories: TFs, water, vegetation, and urban (Fig. 2). To get training and validation data, random stratified sampling was used. A total of 2860 samples were utilized, comprising 65 samples for each feature from 11 different times of images.
This study presents a hybrid framework for TF mapping that integrates the advantages of MSIC, Quantile Composite (10th percentile), and Otsu-based land-water segmentation. In contrast to other methodologies that utilize a singular tidal condition or a single index (e.g., Jia et al., 2021; Zhao et al., 2020), our methodology encompasses both extremes of tidal inundation by employing long-term satellite imagery, independent of tidal elevation models. This composite method is scalable and suitable for diverse tidal regimes, consisting of three key steps: (1) Landsat and Sentinel long-term imagery were obtained and processed using the GEE cloud platform. This procedure comprised first screening for cloud cover, then using QA bands to remove clouds. (2) Then, the time-series change patterns of TFs were analyzed by computing optimal spectral indices (e.g., modified normalized difference water index [mNDWI], normalized difference water index [NDWI]) and applying the maximum value composite algorithm (MVCAmax) to identify maximum water extent and the minimum quantile composite algorithm (MQCAmin) to detect minimum water levels. (3) Subsequently, the Otsu automated thresholding algorithm was employed for precise land-water segmentation, followed by a clipping operation to isolate TF areas by subtracting the minimum water extent from the maximum water extent. (4) We refined the TF maps by removing inland features, islands, and aquaculture ponds using a regional connectivity function and mathematical morphology. (5) For the extracted TFs, a qualitative and quantitative accuracy assessment was conducted.
This study defines TFs as intertidal regions that consist of mudflats and sandflats situated between the highest and lowest recorded tidal extents, including areas that are consistently inundated during high tide and revealed during low tide. These regions may encompass sandy, muddy, and rocky substrates, as observed along the beaches of Karachi and Yangma Island (Murray et al., 2019c; Wang et al., 2020). After getting the high tide and low tide images of Landsat series, we did the spectral analysis of different land covers for quantifying the characteristics of TFs. The average spectral reflectance curves for water, TF, vegetation, and urban classes were derived using eleven high-resolution satellite images acquired at different time periods. The number of samples for each class is 715, 710, 628, and 807 respectively (Fig. 3).
Average spectral. Curves of four land cover types by using the lowest tide images of the Landsat series. The term ‘tidal flat’ refers to a fully exposed TF pixel; ‘water’ refers to both inland and seawater; ‘vegetation’ refers to both inland and intertidal vegetation; and ‘urban’ refers to all inland buildings. TF, tidal flat.
Three widely used spectral indices were evaluated for Landsat and Sentinel-2: the mNDWI, the NDWI, and the Normalized Difference Vegetation Index (NDVI), as shown in Table 1. Although mNDWI and NDWI are commonly used in water extraction, this study applies them within a dual-composite framework designed specifically for TF mapping.
Formulas of the spectral indices used in this study
mNDWI, modified normalized difference water index; NDVI, normalized difference vegetation index; NDWI, normalized difference water index.
We performed an index signature study to assess the suitability of these indices for TF extraction. The 220 sample points for each land cover were utilized to examine the reflectance characteristics of each index for composite images for high and low tide to determine the ideal compositing approach for each tidal scenario.
The mNDWI is employed in the MSIC within MVCAmax approach to capture the maximum water extent, reflecting the highest observed tidal conditions. The MSIC selects the maximum mNDWI value for each pixel over 3 years, effectively delineating the highest observed tidal inundation. A rolling 3-year window centered on each target year was employed to ensure adequate tidal variation to deal with cloud-related gaps while preserving temporal consistency. The following formula for the highest observed tide is used:
S1, S2, …, Sn refer to the various spectral indices in the RS satellite image, such as mNDWI and NDVI.
After preprocessing of satellite data, we produce MSIC images based on mNDWI for high water extent (Fig. S2A in Supplementary Materials). Following Zhang et al. (2019)’s hydrological framework, the maximum water extent composite (MVCAmax) derived from mNDWI effectively incorporates both permanent water bodies and intertidal zones while suppressing vegetation signatures (Fig. 4A). This confirms that mNDWI-based MSIC shows clear boundary separation between TFs, water, and vegetated zones.
Histograms of different land covers for high and low water extent. Histograms representing various land cover classes under high-tide (A and B) and low-tide (C and D) composite conditions.
To represent minimum water extent, we applied a 10th percentile NDWI composite, which captures the shallowest water levels and most exposed intertidal zones over the same 3-year period. The following formula is used:
S1, S2, …, Sn refer to the various spectral indices in the RS satellite image, such as NDWI.
NDWI was chosen over NDVI and mNDWI for this phase because of its superior sensitivity to sediment-rich, shallow waters commonly observed during low tide. NDVI was utilized only for comparative analysis and excluded from the final composite generation because of its restricted discrimination ability under low-tide conditions. Figs 4B–4D illustrate that NDVI and mNDWI quantile composites demonstrated considerable overlap between vegetation and TFs, but NDWI provided a clearer distinction. The ultimate minimum water composite (Fig. S2B in the Supplementary Materials) illustrates the least water coverage, omitting vegetated areas and permanent water bodies.
The dual use of mNDWI for high tide (MSIC) and NDWI for low tide (quantile composite) enables a robust depiction of the TF. This method leverages multi-year time series and extracts spectral extremes independent of tidal elevation data, hence facilitating generalizability to data-sparse regions and improving the reproducibility of TF classification.
This study utilized the OA on the MVCAmax and MQCAmin images to generate a binary classification outcome accurately representing the maximum and minimum water extent. The approach systematically determines the threshold distinguishing water from land categories by analyzing imagery histograms to differentiate between water and land regions (Clement et al., 2018; Donchyts et al., 2016). The formula applied is as follows:
Here, ω0 (t) and ω1 (t) are the pixel class probabilities, and μ0(t) and μ1(t) denote the mean intensities of each class. The approach presumes a bimodal distribution, which was validated through histogram analysis (Figs S3A and S3B in Supplementary Materials). In most cases, Otsu successfully distinguished between water and land classes. For a limited number of composites influenced by cloud or turbidity interference, slight human adjustments within a specified range (0.00–0.45) were implemented to ensure uniformity. Following the threshold application on both composite images, we proceeded with the clipping for the preliminary TF, as shown in the flow chart (Fig. 5).
Flow chart of TF extraction.
The extraction findings were refined to reduce noise and minimize fluctuations, thus assuring a more continuous and precise delineation of the TF line. The smoothing process involved the application of a majority filter utilizing morphological operations, alongside the elimination of all clusters of inland water and island pixels that fell below a specified size threshold. The findings presented in Fig. S4 in the Supplementary Materials demonstrate a clear enhancement in delineation following post-processing, which contributes to more precise area estimation and reliable temporal analysis of TF boundaries.
To evaluate the TF maps generated by our proposed approach, we conducted both quantitative and qualitative accuracy assessments. A stratified random sampling strategy was utilized for the quantitative assessment, leveraging GEE to ensure a balanced representation of TF and non-tidal flat (non-TF) classes. A total of 459 validation points were generated across both study areas. These points were then visually interpreted using low-tide Google Earth imagery to determine their ground truth classification. A confusion matrix was created, and standard classification metrics were computed, including user accuracy (UA), producer accuracy (PA), overall accuracy (OA), and the Kappa coefficient, which measures the agreement between observed and predicted categorization. The metrics offer a comprehensive assessment of classification reliability and the dynamics of error propagation (Bargiel, 2017). The binomial sampling formula (Stehman & Foody, 2019) justified that the sample size is:
The reliability of classification over time, 2014–2016, was reinforced through consistency checks with historical global TF datasets (Murray et al., 2019), which served for qualitative validation.
A confusion matrix was constructed (Table 2) using 250 randomly generated points for TFs and 209 points for non-TF within the study area. The results demonstrate a high overall accuracy (OA) of 94.3% and a kappa coefficient of 0.88 with a 95% confidence interval (CI) margin of error of ±4.67%, indicating strong reliability of the TF maps. The UA and PA for TFs were 93.6% and 96.0%, respectively, further validating the classification method’s resilience.
Confusion matrix of TF validation
| Class | TF | Non-TF | Total | UA | |
|---|---|---|---|---|---|
| Map pixels | TF | 234 | 16 | 250 | 93.6% |
| Non-TF | 10 | 199 | 209 | 95.2% | |
| Total | 244 | 215 | 459 | OA = 94.3% | |
| PA | 96.0% | 92.6% | Kappa = 0.88 |
OA, overall accuracy; PA, producer accuracy; non-TF, non-tidal flat; TF, tidal flat; UA, user accuracy.
This algorithm delineated the spatial extent of TFs across 11 temporal phases over 40 years (1984–2024) at a spatial resolution of 30 m using Landsat data for both KPT and Yangma Island, Yantai.
The extent of the TFs was determined by calculating the area from the extracted data, taking into account the pixel size, as depicted in Table 3 and Fig. 7.
Spatial analysis and rates of area change (ha · y−1) of TF from 1986 to 2024, for KPT, Karachi, and Yangma Island, Yantai, respectively. (A and C) Spatial analysis (B and D) Rates of area change (ha · y−1).
TF maps in KPT, Pakistan, from 1986 to 2024. KPT, Karachi Port Trust; TF, tidal flat.
Temporal variation of TFs area with CI error (ha)
| Year_Range | KPT_Area | KPT_CI_Error | Year_Range | Yangma_Area | Yangma_CI_Error |
|---|---|---|---|---|---|
| 1984–1986 | – | – | 1984–1986 | 972.04 | 45.39 |
| 1986–1988 | 1678.81 | 78.4 | 1887–1989 | 517.06 | 24.15 |
| 1989–1991 | 1651.07 | 77.1 | 1990–1992 | 450.83 | 21.05 |
| 1992–1994 | 1226.50 | 57.28 | 1993–1995 | 415.44 | 19.4 |
| 1995–1997 | 1146.25 | 53.53 | 1996–1998 | 383.23 | 17.9 |
| 1998–2000 | 1136.97 | 53.1 | 1999–2001 | 308.74 | 14.42 |
| 2001–2003 | 947.26 | 44.24 | 2002–2004 | 308.90 | 14.43 |
| 2004–2007 | – | – | 2005–2007 | 372.53 | 17.4 |
| 2008–2010 | 1338.51 | 62.51 | 2008–2010 | 351.64 | 16.42 |
| 2011–2012 | – | – | 2011–2012 | – | – |
| 2013–2015 | 1344.84 | 62.8 | 2013–2015 | 298.04 | 13.92 |
| 2016–2018 | 1396.80 | 65.23 | 2016–2018 | 299.01 | 13.96 |
| 2019–2021 | 1394.76 | 65.14 | 2019–2021 | 258.25 | 12.06 |
| 2022–2024 | 1408.17 | 65.76 | 2022–2024 | 290.12 | 13.55 |
CI, confidence interval; KPT, Karachi Port Trust; TFs, tidal flats.
Table 3 and Fig. 6 illustrate that the TF area adjacent to KPT, Pakistan, underwent substantial transformations from 1986 to 2024.
This region primarily consists of TFs, mangroves, and muddy zones. The TF in the KPT region was categorized into three distinct phases: minimal reduction (1986–1991), reduction (1991–2003), and expansion (2007–2024). In the minimal reduction phase (1986–1991), the TF area was slightly reduced to 27 ha by 1991, indicating a reduction rate of 1.6%.
In the reduction phase (1991–2003), the TF area diminished by 703.81 ha, as illustrated in Fig. 6A, indicating a reduction rate of 43.06%. The mean annual variation during this period was –59.25 ha, accompanied by an average annual decline rate of 3.59%. The most substantial decrease occurred between 1991 and 1994, as illustrated in Fig. 6B, resulting in a 34.6% reduction rate and a decrease of 424.57 ha in the TF area by 1994.
During the expansion phase (2007–2024), the TF area increased by 460.90 ha, marking a growth rate of 48.75%. The average annual change amounted to 22.3 ha, with an average annual growth rate of 2.37% during this period. A notable increase in fluctuation occurred between 2003 and 2010, with a growth rate of 42.5%, as shown in Fig. 6B, where the area expanded from 940 ha to 1338 ha by 2010.
Furthermore, Fig. 6A demonstrates that in the last 40 years, the total area of TFs in KPT along the east coast has diminished by 270 ha, reflecting a total percentage loss over the entire period of 16.1%. The mean annual change during this timeframe was –7.12 ha, corresponding to an average annual reduction rate of 0.41%. The observed fluctuations, characterized by a weak trend (R2 = 0.0046), confirm that the changes are episodic and non-linear, dominated by specific periods of rapid loss and gain, rather than a steady, monotonic trend. The ±4.67% CIs (Table 3, Fig. 4) indicate classification uncertainty, thereby reinforcing the reliability of trend interpretation. The intertidal zone area is progressively increasing toward the southern side of the study region, with the expansion directed seaward.
The extraction of information acquired from RS and the area fluctuations of the tidal zone of Yangma Island, Yantai, are presented in Table 3 and Fig. 8 for the period from 1984 to 2024.We classified the changes in the TF area into two separate periods by meticulous calculations and statistical studies: A reduction phase (1984–1989) and a consecutive reduction phase (1990–2024). These phases highlight the persistent decline in the intertidal zone area over the study period. Yantai’s intertidal zone faces ecological pressures from urbanization, marine aquaculture, and industrial development.
TF maps in Yangma Island, Yantai, China, from 1986 to 2024. TF, tidal flat.
As a popular tourism destination, the development of tourism facilities has disrupted the beach’s natural cross-sectional balance and intensified erosion, highlighting how tourism growth is adversely impacting the intertidal zone. In the initial phase (1984–1989), the TF area experienced a significant decrease of 455 ha, representing a reduction rate of 46.81% (Fig. 6C). This marked a notable decline as shown by Landsat_TM imagery (Fig. S5 in the Supplementary Materials), characterized by a fluctuating trend from 1984 to 1989, with an average annual loss of approximately 28.43 ha (Fig. 6D), corresponding to an average annual reduction rate of 7.8%. TFs were primarily situated on the southern side of the island during this phase, as illustrated in the accompanying map, Fig. 8. Following the reduction phase, the TF area experienced a continued decline of 270 ha, which represents a reduction rate of 43.90% from 1989 to 2024 (Fig. 6C). The average annual change over this period was -6.67 ha (Fig. 6D), yielding an average annual reduction rate of 1.3%. Notably, the decrease in the TF area was predominantly observed in Fig. 8 on the eastern side of the island from 2008 to 2024. In contrast, the southern side initially underwent a reduction until 2010, after which the area began to show a gradual recovery.
Throughout the 14th Five-Year Plan era, spanning from 2021 to 2025, Yantai’s strategic coastal position prioritizes the marine economy and the protection of coastal ecosystems. According to the Yantai Ecological Bureau, in 2020 strictly protected natural coastlines by implementing the Yantai Coastal Spatial Utilization and Protection Plan, establishing a zoning system, and clarifying coastline functions.
To evaluate the temporal persistence and spatial variability of TFs, we performed a frequency-based analysis utilizing pixel-wise aggregation of binary classifications of TFs across all composited time intervals. Each study site acquired a frequency value for every pixel, based on the number of occurrences. It was categorized as a TF throughout the observation period, 11 composited intervals for KPT and 14 for Yangma Island. A frequency value of 1 signifies infrequent occurrence whereas elevated values suggest persistent TF existence across several years. The frequency maps are depicted in Figs 9A and 9B. In the KPT region, 475.5 ha displayed a frequency of 1, signifying the presence of transient or marginal TFs. In contrast, pixels with elevated frequency values, specifically those ≥8, were predominantly located in the southern regions of the research area, indicating morphological stability and recurrent exposure at low tides.
Frequency map of TFs for 11 composite raster’s from 1986 to 2024. (A and B) Frequency map of KPT and Yangma Island, respectively. (C) Sentinel-2 images (2024) for major sites of frequency change in TFs for KPT. (D) Frequency map of TFs for the major sites.
The frequency map delineates notable morphological changes at two key sites: the first is Chinna Creek, located near the coastline, while the second lies east of the Karachi Harbor Entrance (Fig. 9C). These changes are mainly due to human activities, such as port construction, land reclamation, and also natural forces, such as monsoon winds and tidal currents. The study (Fig. 7) illustrates that both areas have undergone processes of accretion and erosion over time. Chinna Creek expanded between 1986 and 2000, but after 2000, it began to shrink, primarily from the north to the south, due to land reclamation activities. The frequency distribution map (Fig. 9D) demonstrates reduced TF visibility, suggesting pronounced anthropogenic influence on the morphological dynamics of these intertidal zones.
The analysis also shows a notable increase in accretion in the East Karachi Harbor between 2008 and 2024. The frequency map (Fig. 9D) indicates a low frequency of TFs in this area, which suggests that TF visibility was limited and likely influenced by significant external factors. During this period, five TF composite images were analyzed, resulting in a frequency scale. The intertidal zone at a frequency of 1 covers an area of 120.16 ha, but as the frequency increases from 5 to 11, this area drastically decreases, approaching nearly zero.
On Yangma Island, most TF pixels had low to moderate frequency values, indicating dynamic exposure patterns influenced by natural forces and coastal development. A compact yet persistent core of 41.36 ha exhibited the highest frequency (F = 14), primarily situated on the southern periphery of the island. These regions are presumably less affected by urbanization. Conversely, lower frequency values on the eastern side correlate with areas impacted by tourism development and marine aquaculture, where TFs are periodically inundated or damaged.
For a detailed analysis of TF dynamics in our study area, we utilized Sentinel-2 satellite imagery spanning the period from 2017 to 2024. Table 4 and Fig. 10 demonstrate that the TF area in the KPT region experienced a slight increase from 2017 to 2024. The low value of R2 suggests a small portion of variability impacted by both stochastic human factors and interannual tidal oscillations, which inherently result in annual variations in TF extent. Conversely, the Yangma Island region in Yantai underwent a decline during the same timeframe.
Area and trend line analysis by Sentinel-2 data for 2017–2024. (A) For KPT. (B) For Yangma Island, Yantai.
Annual TF area with 95% CIs (2017–2024)
| Year | KPT area | KPT CI error | Yangma area | Yangma CI error |
|---|---|---|---|---|
| 2017 | 1492.23 | 69.69 | 291.93 | 13.63 |
| 2018 | 1624.92 | 75.88 | 299.48 | 13.99 |
| 2019 | 1605.06 | 74.96 | 401.01 | 18.73 |
| 2020 | 1227.01 | 57.3 | 438.34 | 20.47 |
| 2021 | 1489.75 | 69.57 | 232.22 | 10.84 |
| 2022 | 1498.94 | 70.0 | 351.32 | 16.41 |
| 2023 | 1236.24 | 57.73 | 278.8 | 13.02 |
| 2024 | 1492.28 | 69.69 | 290.91 | 13.59 |
CIs, confidence intervals; KPT, Karachi Port Trust; TF, tidal flat.
The annual trend analysis for KPT reveals significant fluctuations in TF areas from 2017 to 2024, Fig. 10A. From 2017 to 2018, the area expanded by 132.69 ha, indicating a growth rate of 8.89%. From 2018 to 2020, a significant decrease was recorded, especially between 2019 and 2020, where the reduction rate reached 23.55%, followed by a sharp recovery in 2021 (+21.41%). Following a temporary plateau, there was a second decline in 2023 (–17.53%), followed by another increase in 2024 (+20.71%).
The intertidal zone of Yangma Island exhibited a decreasing trend with a negative slope (R2 = 0.0402) from 2018 to 2024, as illustrated in Fig. 10B. From 2017 to 2024, the TF region of Yangma Island exhibited significant fluctuations. The area increased from 291.93 ha in 2017 to 438.34 ha in 2020, with the most significant growth (33.90%) occurring between 2018 and 2019. Subsequently, there was a significant decrease to 232.22 ha in 2021, representing a 47.02% loss, potentially attributable to land reclamation or tidal conditions. A partial recovery transpired in 2022, attaining 351.32 ha; however, this was not maintained, as the area dropped once more to 278.80 ha in 2023, before seeing a minor increase to 290.91 ha in 2024.
Comparative analysis reveals that Sentinel-2 MSI imagery (10-m) resolves finer tidal creek morphologies than Landsat-8 (30-m), with visual validation confirming accurate detection of >85% of channel features down to 15-m widths. This enhanced capability stems from Sentinel-2’s superior spatial resolution and spectral characteristics. In 2017, the TF area assessed from Landsat-8 was 1267 ha and 275 ha, while Sentinel-2 yielded a larger estimate of 1492 ha and 470 ha for the KPT and Yangma Island, respectively.
To validate the results, we overlaid the Sentinel-2 MSI-derived intertidal area on a single lowest tide image from 2017 for study regions (Fig. 11). Sentinel-2 effectively captured the exposed TFs and minor creeks, demonstrating its capability to accurately map intertidal zones.
Sentinel-2 satellite imagery of the study areas. (A) TF for 2017 in KPT. (B) TF for 2017 in Yangma Island, Yantai. KPT, Karachi Port Trust; TF, tidal flat.
With automated processing on GEE, this study offers a reliable pixel-based compositing approach for multi-decadal TF mapping, producing outputs with constant 30-m (Landsat TM/OLI) and 10-m (Sentinel-2) resolution. High-resolution TF chronologies covering the years 1984–2024 can be reliably generated using the methodology, which is especially useful for capturing fine-scale intertidal dynamics in a variety of coastal environments.
First, the open-access nature and free availability of Landsat and Sentinel-2 imagery have facilitated high-temporal-frequency, high-quality observations, enhancing the ability to monitor the dynamic intertidal zone. Extended observation periods provide a comprehensive dataset that facilitates the analysis and comparison of temporal stages. Observational frequency is critical for generating reliable TF maps, as these efforts rely on the availability of imagery captured near the peak tidal phases. This temporal compositing method improved detection during extreme tidal events, which is essential for accurately delineating TFs in the absence of tide-level records, improving long-term environmental assessment in dynamic coastal regions.
By utilizing the MVCAmax and MQCAmin methods in conjunction with high-quality pixel composites and quantile synthesis, this approach provides an efficient and robust method for identifying imagery corresponding to the highest and lowest tides, respectively. This method, enhanced by the Otsu algorithm, was validated through histogram analysis, indicating bimodal distributions (Figs S3A and S3B in Supplementary Materials). Following thresholding, morphological filtering, and connectivity analysis were applied, which improve boundary continuity and reduce overestimation.
Uncertainty emerged due to lingering cloud cover, discrepancies in tidal stages within fixed composites, inconsistencies among sensors (TM, ETM+, OLI), and spectral confusion involving mudflats or shallow water. The results indicated an overall accuracy of 94.3% and a kappa coefficient of 0.88, accompanied by a 95% CI margin of error ±4.67%, thereby affirming the reliability of the mapping approach.
First, with long-term multi-temporal satellite imagery, capturing images at the peak of high and low tide is difficult, especially in complex tidal regimes. In microtidal habitats like the Yantai coast, low tidal amplitude minimizes spectral difference between wet TFs and neighboring land or aquaculture regions. Second, constant cloud shadows and thin haze, especially during monsoon seasons over Karachi, can contribute to misclassification noise and affect spectral clarity in composite images despite extensive cloud masking utilizing Landsat and Sentinel-2 QA bands. Third, highly turbid waters in deltaic zones, such as Chinna Creek, can replicate the spectral characteristics of wet TFs, resulting in overestimations during specific periods. The observed error margins (±4.67%) in area estimates, as shown in Table 3 and Fig. 6, result from these factors and must be taken into account when analyzing long-term change trends at both study sites.
To validate the accuracy and reliability of our TF maps, we conducted a comparison with existing global-scale TF datasets, specifically the University of Queensland Dataset (UQD). For China, numerous studies have been conducted to map TFs, such as those by Wang at Fudan University and the University of Oklahoma (FUDAN/OU), applied decision tree algorithms on pixel-based and phenology-based approaches to generate annual TF maps at a 30-m spatial resolution for the period 1986–2016. However, there is currently no specific TF mapping available for Pakistan’s coast at the regional or national scale, making the global UQD dataset the only resource for intertidal zone coverage in this region. We evaluated the precision and scope of our maps by comparing the TF areas produced by our methodology with the UQD TF maps from 2008 to 2010, as illustrated in Fig. 12. Our TF area estimates for the KPT and Yangma Island are 1338 ha and 351.6 ha, respectively, by using the Landsat imagery for 2010.
Comparative analysis. (A) KPT, Pakistan, and (B) Yangma Island, China. For KPT, the comparison utilizes the UQD dataset (Murray et al., 2019) for 1989 and 2010. For Yangma Island, the comparison is between the UQD dataset for 1989 and the WF_TF dataset (Wang et al., 2020) for 2016. The red-highlighted area is the most noticeable area of change. KPT, Karachi Port Trust; TF, tidal flat; UQD, University of Queensland Dataset.
In contrast, the UQD reports 1661 ha and 808.2 ha using the Landsat imagery for the same regions and duration. The substantial differences in these area estimates can be largely attributed to the distinct methodologies used in determining tide-level information. The UQD dataset employs a random forest algorithm with a large number of training samples to generate 30-m resolution TF maps within a machine-learning framework. However, this approach frequently misclassifies aquaculture ponds, inland lakes, and salt marshes as TFs due to their overlapping spectral and temporal signatures, leading to overestimation in some areas.
Our methodology employs a dual-composite strategy designed to identify the maximum and minimum tidal extents independently of external tidal models. The highest tide image was generated where the pixel values from multiple scenes were selected to represent extreme tidal conditions based on the mNDWI, which was chosen for its sensitivity to water features. For the lowest tide image, we applied the quantile method (10% quantile), which performs well in separating TFs from vegetation. This method minimizes the misclassification of non-tidal features, such as inland reservoirs and ponds. The method for the minimal water extent ensures that the lowest water levels are represented accurately, even if they fall outside the typical tidal range. It helps in identifying tidal creeks, which improves the overall accuracy of the TF map. However, using higher quantiles, such as WF_TF (Wang et al., 2020) uses the 95% quantile for the highest tide image, could fail to capture the full extent of tidal water coverage, as these values often fall outside the maximal water extent. Therefore, we opted not to use the quantile method for the highest tide image.
In coastal areas, human activities are crucial in determining the extent, erosion and growth of TFs and their impact on TFs is both substantial and inevitable. In our study, spatial analyses and frequency mapping (Fig. 9) reveal distinct morphological changes at both the KPT and Yangma Island locations influenced by accretion, erosion, and human alterations.
TFs in China Creek varied in continuity from 1986 to 2024. In Fig. S6A in Supplementary Materials, almost 100 ha of TFs occurred only once, indicating occasional or unstable exposure, whereas only 30.9 ha lasted throughout all 11 composite periods. This suggests a net decline in long-term stability, with 82% of TF area confined to frequencies 1–5, implying episodic exposure from tidal regime variability and anthropogenic activities.
The TF dynamics on the eastern side of Karachi Harbor have been significantly affected by extensive human activities from 2008 to 2024. In Fig. S6B in the Supplementary Materials illustrates that 120.9 ha of TFs were observed, representing 70% of the total intertidal zone identified during this period. TFs with a frequency of 5, which were consistently visible across all compositing windows, covered <2 ha. The significant reduction in TF recurrence aligns with recorded construction activities, such as the initiation of the Pakistan Deep Water Container Port (PDWCP) in 2009 (Alamgir et al., 2023), along with subsequent dredging and land reclamation efforts. The advancements intended to facilitate large-scale commercial shipping have resulted in habitat disruption, alterations to coastlines, and a reduction in the persistence of TFs. Karachi Coastal Comprehensive Development Zone (KCCDZ) and Sindh Coastal Development Authority (SCDA) manage development projects in the coastal regions of Karachi (Khan et al., 2024). It is involved in monitoring coastal development and promoting sustainable practices along the coastline.
On Yangma Island, years of coastal industrialization, aquaculture, and tourism development have resulted in significant TF degradation. Since 1987, the region has undergone significant fragmentation of TFs, especially on the island’s eastern front. Tourism-related development has also disturbed natural coastline dynamics, increasing susceptibility to erosion. Nevertheless, Protective measures established under the 14th Five-Year Plan (2021–2025), such as zoning schemes and the Blue Bay restoration initiative, are starting to alleviate degradation and restore ecological balance in specific areas (Hepburn et al., 2021).
The emergence of advanced RS has markedly improved our capacity to observe and assess coastal ecosystems, particularly TFs, which are highly dynamic and subject to frequent environmental fluctuations. Recent developments in machine learning and automated processing algorithms have optimized the extraction of TF features from RS data, improving accuracy and efficiency, which is affected by the tidal dynamics (Avtar et al., 2020). Our pixel-based composite imagery method can enable detailed estimation of tidal extent utilizing high spatial resolution satellite imagery, such as Sentinel-2, with microwave data from Sentinel-1 (Veloso et al., 2017), as well as MODIS and WorldView-3 data, overcoming the difficulties imposed by cloud cover is necessary for future research. This method can also be applied to different wetland features by applying different spectral indices. By incorporating diverse data types, such as tide measurements from the tidal stations, high-resolution topographic models, and environmental variables, these integrated approaches can provide a more comprehensive view of TF dynamics. A promising future direction is to scale the proposed pixel-based composite approach to generate a national-scale TF inventory for Pakistan and China, providing a consistent and accurate baseline for macro-level coastal management.
Previous studies on TF mapping have primarily relied on supervised or unsupervised classification methods, composite techniques, or auxiliary data, which often face limitations in spatial and temporal resolution. This study addresses a rapid and automated approach to extracting precise TF data and the influence of human activities on TFs.
- 1)
By leveraging the cloud computing capabilities of the GEE platform and utilizing a pixel-based composite algorithm, we generated high-resolution TF maps (30 m and 10 m) for KPT, Pakistan, and Yantai, China, using Landsat TM/OLI imagery for a 3-year composite and Sentinel-2 data for annual mapping. Our methodology achieved high accuracies (overall accuracy = 94.3%, Kappa = 0.88, and CI = ± 4.67%), outperforming conventional datasets.
- 2)
The study utilized the mNDWI combined with the MVCAmax to identify maximum water extent and the NDWI with the MQCAmin to determine minimum water extent. For enhanced precision, the Otsu algorithm and post-processing techniques were applied to delineate the spatial extent of TFs.
- 3)
Over 40 years, the 30 m spatial resolution maps revealed three distinct phases of TF dynamics along the KPT Karachi: minimal reduction (1986–1991), reduction (1991–2003), and expansion (2003–2024). During the past 4 decades, TF areas in the KPT region decreased by 270 ha, reflecting a reduction rate of 16.1% with an average annual change of –6.75 ha. In contrast, for Yangma Island, Yantai, the TF area underwent a drastic reduction over two phases: the first reduction phase (1984–1989) was followed by a continued decline from 1990 to 2024, with a total decrease of 682 ha, representing a 70.16% reduction. Annual Sentinel-2 monitoring from 2017 to 2024 acquired precise area variations and validated spatial accuracy through visual assessment.
This research offers an achievable approach for the long-term monitoring of TFs in coastal environments. The method’s independence from tidal models and training samples reduces overestimation common in supervised classifications, providing a more accurate baseline for management. The resulting maps offer reliable data for TF management, sustainable coastal zone development, and supporting related sustainable development goals (SDGs) and scientific research.