Exploring water accumulation dynamics in the Pearl River estuary from a Lagrangian perspective

Investigating water accumulation is essential for understanding estuarine mass distribution and ecosystem management. In this study, we examined the water accumulation dynamics in the Pearl River estuary (PRE) from a Lagrangian perspective. Generally, there is a notable negative correlation between the horizontal velocity divergence (∇_h⋅V_h) and the accumulation. Influenced by density fronts and velocity convergence, significant bottom-layer accumulation of particles was observed in the western part of the estuary and the Hong Kong waters during summer, whereas this accumulation moved landward in winter. Subregions with distinct accumulation patterns and interconnections were identified and combined with the trajectories. In summer, the western part of the estuary and the Macau waters exhibit significant net negative ∇_h⋅V_h values, and strong density fronts are major accumulation targets, attracting particles from the entire estuary. Conversely, the eastern part of the estuary and Hong Kong waters exhibit significant westward motion, influencing the western side. In winter, particles are more likely to accumulate in their initial locations. The upper part of the estuary becomes a major accumulation area due to the obstructive density front and decreased river discharge. Tidal currents and river discharge mainly control water accumulation in the estuary by changing the mixing or current intensity. The weakening of tidal currents and river discharge induce intensified bottom intrusion and the landward movement of accumulation.

1 Introduction

The Pearl River estuary (PRE), located in the northern South China Sea (NSCS) (Fig. 1a), is influenced by the East Asian monsoon, with northeasterly winds prevailing in winter and southwesterly winds prevailing in summer (Li and Li, 2018). Thus, in the PRE, winter is characterized as a dry season, while summer is characterized as a wet season due to the heavy rainfall induced by the moist air from the South China Sea. Consequently, river discharge in summer (∼20 000 m³ s⁻¹) is approximately 5 times greater than that in winter (∼3600 m³ s⁻¹) (Harrison et al., 2008). This is in contrast to many other river deltas, such as the Mississippi Delta, where river discharge reaches a maximum in winter and spring but decreases in summer and fall (Lane et al., 2007).

As the PRE is a bell-shaped estuary, its width increases from approximately 5 km at the upper end to 35 km at the lower end. Despite having two narrow, deeper channels (∼20 m in depth), the PRE is shallow, with a water depth of approximately 2–10 m. The PRE is a partially mixed estuary in which circulation is jointly controlled by river discharge, tides, wind, and topography (Ascione Kenov et al., 2012; Banas and Hickey, 2005; Gong et al., 2009; He et al., 2023, 2022; Liu et al., 2020). There are two distinct dynamic regimes in the PRE. The narrow upper part of the PRE exhibits classical gravitational circulation, while in the wider lower part of the PRE, where the Coriolis effect becomes significant, the interaction between the monsoon-driven shelf current and the topography complicates the circulation (Dong et al., 2004; Wong et al., 2003; Zu and Gan, 2015). Gravitational circulation occurs in the two deep-channel regions, while currents exhibit precise seasonal characteristics across the shallower western part of the estuary. Geostrophic wind-driven coastal currents intrude into the PRE during the summer upwelling season (Zu and Gan, 2015), whereas seaward buoyancy-driven coastal currents flow out of the PRE during winter (Dong et al., 2004; Wong et al., 2003). The alternating spring–neap tide cycle and variations in river discharge play crucial roles in modulating stratification and mixing within the PRE (Mao et al., 2004; Pan et al., 2020; Zu et al., 2014). Strong tidal mixing in the middle of the PRE has led to the conversion of estuarine river plumes into buoyancy-driven coastal currents (Dong et al., 2004; Zu et al., 2014).

In recent decades, the state and sustainability of estuaries have been increasingly threatened by environmental issues caused by anthropogenic activities and climate variability (Dai et al., 2006; Fok and Cheung, 2015; Shen et al., 2023). Estuarine circulation and the associated mass transport are critical for addressing these environmental challenges. Specific areas of water accumulation and the extent of their exchange with the open sea significantly impact the dispersion of oceanic pollutants and biogeochemical conditions (Acha et al., 2004; Hinojosa et al., 2011). In estuaries, some regions are more likely to attract water due to complicated current circulation, making them stronger horizontal-convergence targets for certain materials (Wang et al., 2022). For example, the salt wedge acts as a significant pollutant sink in an estuary, with the salt wedge in the Río de la Plata estuary consistently exhibiting higher concentrations of microplastics (Acha et al., 2003; Vermeiren et al., 2016). Areas with higher concentrations of nitrogen and phosphorus are usually eutrophic (Tao et al., 2021). Heavy-metal pollution in estuaries has been observed in areas with concentrations of fine particles (Balachandran et al., 2005). Therefore, identifying accumulation areas in estuary–shelf systems is crucial for accurately surveying pollutant sinks (Mestres et al., 2006; Tao et al., 2021; Vermeiren et al., 2016; Wang et al., 2016).

With intensified human activities, pollutant sinks related to accumulation phenomena in the PRE have garnered attention. Tao et al. (2021) revealed that the upper part of the PRE is a key sink for nitrogen and silicate. Zhang et al. (2013) found that trace elements tend to accumulate on the PRE's western side. Higher concentrations of microplastics have been observed in the western part of the estuary and the Hong Kong waters (Lam et al., 2020). Similarly, studies on hypoxia have shown that the convergence of buoyancy-driven currents and wind-driven shelf flows contributes to the formation of stable water columns, creating favorable conditions for the development of hypoxic zones (e.g., Li et al., 2021, 2020). The high frequency of hypoxia in the estuary during summer is linked to the strong stratification of the water column (Cui et al., 2019; Zhang and Li, 2010). These studies on accumulation patterns in the PRE focus on measuring pollutant concentrations (Tao et al., 2021), estimating pollutant accumulation rates (Zhang et al., 2009), and discussing pollutant sources (Ye et al., 2012), but they lack a comprehensive understanding of spatial accumulation patterns and the underlying physical controls.

The Lagrangian tracking method provides an effective method for analyzing the transport processes in physical oceanography (Jalón-Rojas et al., 2019; Van Sebille et al., 2018). It involves analyzing the sizable virtual particle trajectories calculated from simulated three-dimensional Eulerian and time-varying velocity fields, which capture complex real-world dynamical processes. Its applications extend across various domains, including interocean exchange (Haza et al., 2016), pathway analysis (Jalón-Rojas et al., 2019), and the impact of ocean currents on ecosystems (Chenillat et al., 2015; Dawson et al., 2005; Lebreton et al., 2012; Paris et al., 2012).

In this study, Lagrangian tracking and analysis were utilized to examine the features of the accumulation regions in the estuary–shelf system of the PRE and to explore the role of various forcing factors. The results enhance our understanding of the environmental effects of multiscale processes. The remainder of this paper is organized as follows. Section 2 introduces the Lagrangian model and numerical solutions of the Markov chain. Findings regarding the accumulation region and transport connections in the PRE are discussed in Sect. 3. Section 4 presents the roles of hydrodynamic factors. Finally, Sect. 5 presents the conclusions of this study.

Figure 1 (a) The location and topography of the Pearl River estuary (PRE) and the adjacent shelf. (b) The bathymetry of the PRE and the adjacent shelf. The dotted and dash-dotted black lines represent the 25 and 50 m isobaths, respectively. The yellow line defines the seaside boundary of the PRE. LDY, MO, HK, and NSCS represent Lingdingyang, Macau SAR, Hong Kong SAR, and the northern South China Sea, respectively. The red line represents the model boundary.

2 Methods

2.1 Numerical simulation and Lagrangian particle tracking

In this study, multiscale circulation in the PRE and adjacent shelf regions was simulated using the Regional Ocean Modeling System (Shchepetkin and Mcwilliams, 2005). The model region covers the estuary and adjacent shelf regions between 112.3–115.68° E and 20.89–23.13° N (Fig. 1b). An orthogonal curvilinear grid was employed, with the resolution increasing gradually from approximately 1 km over the shelf to ∼200 m within the estuary. In the vertical direction, we used the terrain-following s coordinate (Song and Haidvogel, 1994) to discretize the water column into 30 layers.

Monthly averaged data on riverine discharge were obtained from the Ministry of Water Resources of China. River discharge was approximately 30 000 and 10 000 m³ s⁻¹ in summer and winter, respectively (Fig. S1 in the Supplement). Wind forcing, heat flux, and precipitation data were obtained from ERA5 (European Centre for Medium-Range Weather Forecasts Reanalysis v5) atmospheric-reanalysis data, provided by the European Centre for Medium-Range Weather Forecasts (ECMWF), and were used to force ocean circulation through the implementation of the bulk computation algorithm (Fairall et al., 2003). The shelf current, obtained from a coarser model with strong validation covering the NSCS (Fig. S2), provided information on barotropic and baroclinic velocities, temperature, salinity, and sea levels along the boundaries of the PRE (Deng et al., 2022). The vertical turbulence and diffusion coefficient were determined using the Mellor–Yamada 2.5 turbulence closure module (Mellor and Yamada, 1982), which provides the turbulent-mixing coefficient. Simulations of tidal currents, accounting for forcing along open boundaries, were obtained from Zu et al. (2008), including the tidal harmonic constants (M₂, S₂, K₂, N₂, K₁, O₁, P₁, Q₁, and M₄). This model, based primarily on climatological data, was carefully verified using satellite remote sensing and long-term observations to ensure an accurate representation of the hydrodynamic properties (Figs. S3 and S4). Overall, the model accurately captured the seasonal variability in the hydrodynamic features in this region, and it has been used in previous studies (Cai et al., 2022; Chu et al., 2022b; Cui et al., 2024).

Particle trajectories were traced using the three-dimensional offline Larval TRANSport Lagrangian model (LTRANS v.2b), which captures complicated dynamical processes in the real world through Eulerian flow fields and turbulent mixing from the hydrodynamic model (Chu et al., 2022a; Liang et al., 2021; North et al., 2011). To reasonably calculate Lagrangian trajectories within estuary–shelf circulation systems, the tracking model considers advection, turbulence, individual particle behaviors (e.g., vertical sinking, floating, or swimming velocity), settlement, and boundary behaviors during particle trajectory simulations. The fourth-order Runge–Kutta scheme was implemented to handle the advective terms and yield accurate particle trajectories (Dippner, 2009). Considering the random walk of water parcels in the ocean, the model adopts different diffusivity coefficients to control vertical and horizontal turbulence (North et al., 2006; Zhong and Li, 2006). In the surface-tracking (bottom-tracking) case, 8386 particles were uniformly distributed at the water surface (along the water bottom) across the estuary and the adjacent shelf with a 0.01° interval. Particles were released every 2 d and tracked for 20 d. The hydrodynamic simulation results were stored every 20 min in January and July to drive particle tracking during summer and winter, respectively. During trajectory tracking, the time step for particle tracking was 30 s, and the locations of the particles were recorded every 20 min.

2.2 Markov chain

In this study, we investigate accumulation features and connectivity across different parts of the domain using Markov chain statistical analysis. This approach describes the future state of a random variable based on the particle's initial state, allowing us to compute a proper transition matrix (Drouin and Lozier, 2019; Drouin et al., 2022; Miron et al., 2021, 2019, 2017; Van Sebille et al., 2018). To this end, we divided the study area into rectangular grids with dimensions of 0.1° and then calculated the probability of particles moving between different grids within the time interval of dt as follows:

$$\begin{array}{}\text{(1)}& p^t={\displaystyle \frac{n_{ij}^{t_{\mathrm{0}}+\mathrm{d}t}}{n_i^{t_{\mathrm{0}}}}},\end{array}$$

where $n_i^{t_{\mathrm{0}}}$ represents the number of particles released in grid i at initial time t and $n_{ij}^{t_{\mathrm{0}}+\mathrm{d}t}$ represents the number of particles arriving at grid j from grid i after interval time dt. Thus, p^t varies for each different release grid i, arrival grid j, and period (t₀, t₀+dt). Subsequently, the evolution of the initial distribution ($D^{t_{\mathrm{0}}}$) into a future state ($D^{t_{\mathrm{0}}+\mathrm{d}t}$) is achieved using vector–matrix multiplication as follows:

$$\begin{array}{}\text{(2)}& D^t=D^{t_{\mathrm{0}}}\times \prod p^t.\end{array}$$

In this study, $D^{t_{\mathrm{0}}}$ is defined as a uniform initial mass distribution across the entire domain – i.e., $D^{t_{\mathrm{0}}}=[\mathrm{1},\mathrm{1},\mathrm{1},\mathrm{\dots }\mathrm{1}]$. D^t is the evolution of the initial condition under complicated hydrodynamic motion, which is calculated using the multiplication of $D^{t_{\mathrm{0}}}$ with the transition matrix P^t (Fig. 2). The transition matrix P^t is derived from the multiplication of the probability matrix p^t, illustrating the net probability trend between different regions in the study area. Areas with high D^t values act as strong accumulation targets for particles.

Figure 2 Schematic representation of the calculation of the Markov chains.

The transport matrix generated by Markov chains aids in predicting longer-timescale transport pathways through limited short-time trajectories and is widely used to explore transport pathways and the connectivity of water parcels between different regions. For instance, Van Sebille et al. (2018) investigated the pathways of significant currents in the Agulhas Current, while Drouin et al. (2022) found that there was no effective water transition between the North Brazil Current and the Gulf of Mexico. Studies by Jönsson and Watson (2016) and Drouin et al. (2022) have shown that the choice of dt value does not significantly influence the results of transport patterns. For the present analysis, we selected a time interval (dt) of 2 d.

3 Results

3.1 Accumulation pattern and regional connectivity

We set up a standard real-world case to explore the accumulation pattern and regional connectivity, considering river discharge, tidal currents, and wind forcing. Before examining the detailed transport structure, we examined the offshore transport speed of particles after their release. Incorporating all the particle trajectories, Fig. 3 shows the percentage of particles remaining within the PRE during summer and winter. Here, the seaside boundary of the PRE is defined by a 25 m isobath (yellow line in Fig. 1). If particles move beyond the seaside boundary of the estuary, they may return to the study area due to tidal currents and wind forcing. However, once particles leave the model domain, they will not be backed again. Overall, the decay speed is faster in summer than in winter, with offshore motions at the surface layer always being quicker than those at the bottom layer. Using a value of 20 % as a threshold (dotted green line in Fig. 3), during summer, approximately 80 % of the surface particles exit the estuary in approximately 10 d, whereas the bottom particles take around 15 d. In contrast, during the winter season, with reduced river discharge, the exit time increases to 20 d (for the surface particles) and 25 d (for the bottom particles). Thus, across the entire domain, after being tracked for 20 d, most particles move beyond the seaside boundary, and their trajectories are used in the following analysis. Using the trajectories of the released particles over 20 d, we explored the final evolved state of the particles to quantify the accumulation targets resulting from the complex hydrodynamics of estuarine circulation.

Figure 3 Changes in the percentage of particles remaining in the estuary during tracking.

Examining the circulation in the standard case (Fig. S5), intensified offshore motions are consistently observed at the surface layer in both summer and winter. Due to the reduction in river discharge, the offshore current in the upper part of the estuary is weaker during winter than during summer. The shelf current intrudes into the estuary from the western side of Hong Kong SAR and reaches the middle of the estuary at the bottom, with this bottom intrusion being more pronounced in winter. We select the AB transect (pink line in Fig. 4a) along the intrusion channel to further examine the density stratification and vertical circulation (Fig. 4e and f). A large amount of river discharge in summer pushes fresh water to the seaside, resulting in the seaward movement of seawater (indicated by the isopycnal line marking 1015 kg m⁻³). In winter, reduced river discharge weakens the ability of fresh water to resist the intrusion of saltwater.

Figure 4 illustrates D^t, where regions of high values represent favorable targets for particle accumulation. As a result, compared to particles at the bottom layer, surface particles tend to escape from the estuary quickly due to the strong offshore current at the surface layer, resulting in very low accumulation (Fig. 4a and c). On the contrary, at the bottom layer (Fig. 4b and d), the accumulation targets are much more apparent. In summer, high accumulations are observed on the western and southwestern sides (near Macau SAR) and in Hong Kong waters and are associated with the westward transport current of Hong Kong waters and the intensified intrusion current in the lower part of the estuary. Although these regions are close to the open shelf, they are primary accumulation zones that capture particles that are initially released in other areas of the domain. Cruise observation data (2014–2018) also indicate that these regions generally have a high frequency of hypoxia due to the accumulation effect of hydrodynamic processes (Li et al., 2021).

During winter, with reduced river discharge, gravitational circulation is generally preserved but with weakened intensity. The accumulation regions move shoreward, particularly in the upper part of the estuary, due to decreased circulation and the landward movement of the dense intrusions. Due to intensified northeasterly winds, the westward transport of Hong Kong waters and the adjacent shelf current are strengthened, leading to higher accumulations along the western side of the estuary.

Figure 4 Panels (a–b) show particle mass (D^t; colored areas) at the surface layer and bottom layer, respectively, during summer. Higher values represent stronger accumulation. Panels (c–d) display the same information as panels (a–b) but for winter. Panels (e–f) show the along-transect velocity (colored areas and arrows), with positive values indicating onshore intrusion, and density contours of 1010 and 1015 kg m⁻³ (red lines) in AB during summer and winter, respectively. The location of AB is shown in Fig. 4a.

Compared to that at the bottom layer, the quicker motion at the surface layer does not distinctly reveal an accumulation pattern. Hence, the study focuses on the accumulation pattern and regional connectivity at the bottom layer. Based on the bottom-layer accumulation patterns (Fig. 4b and d), the study area is divided into six subregions: the upper part of the estuary (UPPER), the western part of the estuary (WEST), the eastern part of the estuary (EAST), Hong Kong waters (HK), Macau waters (MO), and shelf waters (SHELF). Subsequently, using the trajectories, the transition matrix (P^t) across each region during the tracking period is examined (Fig. 5). Sensitivity experiments show that different tracking periods do not significantly affect the main patterns observed in the transition matrix.

During summer (Fig. 5a), the initial region consistently contributes significantly to particle accumulation. The WEST and MO regions, which exhibit the most significant attraction ability, accumulate water from almost all subregions, including the shelf. Onshore intrusion from the shelf water, offshore motions from the UPPER region due to river discharge, and westward transport from the EAST and HK regions result in the convergence of water in the WEST and MO regions, facilitating particle accumulation in these regions. In contrast, the EAST and HK regions mainly accumulate particles from the eastern side of the estuary. Despite the Eulerian currents, there is onshore intrusion from the shelf and offshore transport from the UPPER region toward these areas.

During winter, accumulation regions move shoreward and are mainly located in the UPPER, WEST, and MO regions (Fig. 5b). The particles that accumulate in these regions are mainly from the initial-release regions, making it more challenging for water to leave the initial regions. Hong Kong waters can influence almost the entire estuary, especially in terms of contributing to arrival targets in the WEST region.

Figure 5 The connection between six subregions at the bottom layer during summer (a) and winter (b). The x and y axes represent the arrival and release regions, respectively.

3.2 Hydrodynamic control on accumulation

The bottom divergence of the horizontal current, i.e., ${\mathrm{\nabla }}_{\mathrm{h}}\cdot {\mathit{V}}_{\mathrm{h}}=\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}$ (where u and v represent bottom zonal and meridional velocity, respectively), is calculated to examine its influence on the identified bottom-layer accumulation regions (Fig. 6a and b). We established a connection between the average D^t value in each subregion and the divergence of the horizontal current ∇_h⋅V_h. Across the various identified subregions, a significant negative correlation between ∇_h⋅V_h and particle mass is observed, with a correlation coefficient of 0.74 in summer and 0.76 in winter. This suggests that the net negative ∇_h⋅V_h values, i.e., the convergence of water, provides favorable conditions for the accumulation of water and particles. Such patterns are noted in major accumulation regions, such as the WEST, MO, and UPPER regions, during both summer and winter.

It is also noted that the EAST and HK regions, which exhibit particle accumulation, show a weak net positive value of ∇_h⋅V_h. This implies that while ∇_h⋅V_h facilitates accumulation under certain conditions, actual net accumulation should consider the cumulative effects of velocity convergence along the trajectories, which occur at different locations throughout the movement. The spatial distribution of ∇_h⋅V_h also illustrates that intensified negative values occur in regions with high particle mass (Fig. 6c and d), such as the significant negative ∇_h⋅V_h values observed in the MO and WEST regions during summer and in the UPPER region and central part of the HK region during winter. However, the Eulerian perspective of ∇_h⋅V_h presents a complex distribution of alternative positive and negative values, making it challenging to straightforwardly identify areas of net accumulation. In contrast, Lagrangian tracking offers a clearer understanding of these regions because it captures the cumulative effects of water motion over time.

Figure 6 (a–b) Scatterplots illustrating regional D^t against ∇_h⋅V_h for various subregions during summer and winter, respectively. (c–d) Horizontal distribution of ∇_h⋅V_h during summer and winter. ∇_h⋅V_h is normalized using the largest divergence value in this area. Purple lines each represent the isopycnic line of 1010 kg m⁻³.

The salinity front acts as a barrier to particle transport and plays an important role in accumulation regions. For example, coarse particles accumulate at the bottom of the salinity front (Defontaine et al., 2020; He et al., 2018; Vermeiren et al., 2016). The coupling effect of velocity convergence and the front leads to particle accumulation in the middle of the estuary during winter and in the lower part of the estuary during summer (Malli et al., 2022). In summer, heightened river discharge creates a density front, $\left(G=\sqrt{{\left(\frac{\partial \mathit{\rho }}{\partial x}\right)}^{\mathrm{2}}+{\left(\frac{\partial \mathit{\rho }}{\partial y}\right)}^{\mathrm{2}}}\right)$ (Fig. 7), affecting the transport of particles and their accumulation. The location of the front roughly aligns with the outer-boundary region of the profoundly negative ∇_h⋅V_h values across the WEST and MO regions (Fig. 6c), further supporting accumulation in these regions as the front hinders offshore movement. During winter, the intrusion reaches the middle of the estuary, creating a blocking effect on mass transport, which tends to remain confined to the upper part of the estuary (Lima et al., 2015; Lima et al., 2014). This retention likely contributes to the observed accumulation patterns, where particles remain confined upstream due to the blockage effect of the front. ∇_h⋅V_h in winter also reflects this dynamic as regions with negative values move shoreward.

Figure 7 (a–b) The density front (kg m⁻⁴) of the standard case at the surface layer and bottom layer, respectively, during summer. Panels (c–d) show the same information as panels (a–b) but for winter.

4 Discussion

Tides and rivers have an essential influence on the estuarine circulation and associated mass transport in the PRE, whereas winds mainly affect the shelf current and have less influence on mass accumulation within the estuary (figures not shown here). To examine their contributions to the accumulation, we conducted two additional experiments. In the reduced-river-discharge case, the river discharge magnitude in the forcing file is reduced by 20 %. The magnitude of the tidal current is set to zero to remove the tidal currents. The accumulation regions and connections across the different regions are then investigated in the same way using the same Lagrangian tracking and analysis approach.

4.1 Tide

Because the tidal current affects the mixing intensity in the water column, once it is removed, offshore motion in the UPPER region is reduced during both summer and winter. In summer, the westward transport of Hong Kong waters and adjacent shelf currents, flowing south of the WEST region, is strengthened (Fig. S6a). In winter, water from the shelf intrudes directly landward and reaches the upper part of the estuary. Along the AB transect (Fig. 4a), we further examine the changes in the density stratification and vertical structure of circulation (Fig. 8c and d), revealing intensified and landward bottom-water intrusion associated with stronger density stratification. Removing tidal currents leads to weaker mixing and intensified stratification, as observed in previous studies where weaker tidal intensity contributed to increased stratification in the PRE (e.g., Pan et al., 2020).

The D^t anomalies in the tidal-current-removal case and the standard case, i.e., $D_{\text{no-tide}}^t-D_{\mathrm{standard}}^t$, are compared in Fig. 8. Consequently, the removal of tidal forces intensifies particle accumulation in the middle of the WEST region and in the northeastern parts of the EAST region during summer (Fig. 8a). Conversely, in the MO and HK regions, the absence of tides generally reduces the likelihood of accumulation, particularly in the MO region. In winter, intensified onshore intrusion leads to higher accumulation in the UPPER regions (Fig. 8b). In contrast, in the WEST region, tides significantly hinder particles from exiting the estuary due to the intensified landward motion of the westward shelf current, resulting in a negative accumulation anomaly when tides are removed. Thus, the tidal current mainly resists bottom-water intrusion into the midwestern part of the estuary by changing the mixing intensity and density structure.

Figure 8 Panels (a–b) show the particle mass (D^t) anomaly in the tidal-current-removal case during summer and winter, respectively. Negative values represent strengthened offshore transport without tidal currents. Panels (c–d) show the along-transect velocity (colored areas and arrows), with positive values indicating onshore intrusion, and density contours of 1010 and 1015 kg m⁻³ (red lines) in AB during summer and winter, respectively, in the tidal-current-removal case.

Similar to the standard case, there is a negative correlation between the ∇_h⋅V_h values and particle mass D^t (Fig. 9a and b). The arrows represent the changes in ∇_h⋅V_h and particle mass due to the removal of tidal currents in the region, with significant changes in particle mass observed. The decrease (increase) in summer accumulation in the MO region (EAST region) correlates with a weakening (strengthening) of the net convergence of the current due to the absence of tidal currents. Similarly, changes in winter accumulation in the UPPER, WEST, and EAST regions are also related to current divergence. However, it is noteworthy that in the WEST and MO regions, where there is relatively strong ∇_h⋅V_h, changes in accumulation are not significantly influenced by this divergence.

The impact of tides on water transport across the various identified subregions is investigated, as shown in Fig. 9c and d. In summer, removing the tide results in a strengthened attraction ability in the EAST and UPPER regions. The WEST region continues to exhibit a strong converging trend owing to the strengthened landward-intrusion current, attracting particles across the estuary. Notably, in the absence of tides, the MO region no longer attracts particles from other regions. Furthermore, particles originating from the EAST and UPPER regions ceased to reach the HK and MO regions. During winter, the UPPER region shows increased attraction to particles from the WEST, EAST, and MO regions, which is associated with intensified onshore currents from the lower part of the estuary. The converging trend in the WEST region diminishes rapidly, reducing its role as a significant trap for water from other regions. There are no significant changes in the EAST and HK regions; only the MO region traps more water from its own region.

Figure 9 (a–b) Scatterplots illustrating particle mass (D^t) against ∇_h⋅V_h in the tidal-current-removal case for various subregions during summer and winter, respectively. Arrows represent the changes in ∇_h⋅V_h and particle mass due to the removal of tides in each subregion. (c–d) Connections between the six regions at the bottom layer in the tidal-current-removal case during summer and winter, respectively.

Figure 10 Panels (a–b) show the particle mass (D^t) anomaly in the reduced-river-discharge case during summer and winter, respectively. Negative values represent strengthened offshore transport with reduced river discharge. Panels (c–d) show the along-transect velocity (colored areas and arrows), with positive values indicating onshore intrusion, and density contours of 1010 and 1015 kg m⁻³ (red lines) in AB during summer and winter, respectively, in the reduced-river-discharge case.

4.2 River discharge

Based on the circulation pattern resulting from reduced river discharge, a significant decrease in seaward currents is observed throughout the entire estuary during both summer and winter (Fig. S6c and d). Similar to the results shown in Fig. 8, the anomaly results are calculated by subtracting the D^t value of the standard case from the D^t value of the reduced-river-discharge case (Fig. 10a and b). Along the AB transect, landward velocity increases in the upper part of the estuary, while seaward motion accelerates in the lower part of the estuary during winter (Fig. 10d). This change is associated with the divergence and convergence of accumulation in the middle and lower parts of the estuary in winter (Fig. 10b). Increased landward currents result in accumulation in the UPPER region. In summer, reduced stratification leads to a well-mixed water column in the inner part of the estuary, and the landward movement of the density front to the middle of the estuary causes a blockage of water in the upper part of the estuary (Fig. 10c). As a result, accumulation in the inner part of the estuary increases during summer when river discharge is reduced (Fig. 10b).

Generally, major changes in regional accumulation are governed by current convergence. In both seasons, variations in net velocity convergence resulted in corresponding changes in particle mass in the UPPER, WEST, and EAST regions (Fig. 11a and b). Increased particle mass is accompanied by decreased net divergence.

Similarly, Fig. 11c and d illustrate that reduced river discharge greatly impacts the connections between the UPPER, WEST, and MO regions. During summer, in the UPPER region, the converging ability becomes stronger due to the rapid reduction in offshore currents and the blocking of the landward density front. Consequently, more particles from the UPPER region remain in that region and are less likely to move to the EAST and MO regions. Moreover, reduced river discharge results in fewer particles from the MO region being transported into the WEST region. Instead, the EAST region contributes more to the retention of particles in the WEST regions. During winter, the interconnections among the six subregions are similar to those in the standard case (Fig. 5b). It should also be noted that the WEST region experiences lower particle convergence from the EAST and HK regions. Additionally, the MO region attracts particles from the WEST region rather than from the SHELF region.

Figure 11 (a–b) Scatterplots illustrating particle mass (D^t) against ∇_h⋅V_h in the reduced-river-discharge case for various subregions during summer and winter, respectively. Arrows represent the changes in ∇_h⋅V_h and particle mass due to the reduction in discharge in each subregion. (c–d) Connections between the six regions in the reduced-river-discharge case during summer and winter, respectively.

Figure 12 Schematic of the accumulation connections in the PRE during summer (red arrows) and winter (green arrows). The red colors represent the areas of high accumulation in summer, while the green colors represent these areas in winter. The stars indicate where the tide dominates the current, and the triangles represent river discharge.

5 Conclusions

In this study, the Lagrangian method and Markov chains were applied to illustrate the accumulation trends in different PRE regions during typical monsoon seasons.

The accumulation patterns and regional connectivity were derived from the Markov chains. We found that, generally, surface offshore transport is always quicker than bottom offshore transport due to strong offshore currents and river discharge, leading to relatively low levels of accumulation at the surface layer. Conversely, significant accumulation occurs at the bottom, primarily in the lower part of the estuary during summer, and moves shoreward during winter due to reduced river discharge and intensified density front. Based on these bottom accumulation patterns, we identified six subregions in the PRE: UPPER, WEST, EAST, MO, HK, and SHELF. Across these subregions, there is a negative correlation between the net divergence (∇_h⋅V_h) and accumulation, where more intense negative ∇_h⋅V_h values provide favorable conditions for water accumulation. Connections between the six subregions were discussed to illustrate the transport structure in the PRE. During summer, regions like the WEST and MO regions, which are characterized by significant net negative divergence and strong fronts, serve as powerful accumulation targets that attract particles from almost the entire estuary. Conversely, the EAST and HK regions show a westward motion and are transported to the western part of the estuary. In winter, accumulated regions exhibit self-correlations, with particles being more likely to remain in their initial regions. The UPPER region becomes a major accumulation area due to the blocking of the density front and the significantly decreased river discharge. Hong Kong waters are transported throughout the entire estuary, primarily contributing to the accumulation in the WEST region under westward currents.

Sensitivity experiments were conducted to evaluate the impacts of tidal currents and river discharge on accumulation patterns. Generally, tidal currents and river-induced gravitational circulation affect accumulation in different ways, influencing various regions of the estuary. Their joint effects control the accumulation pattern. Tidal currents enhance accumulation in the WEST and UPPER regions during winter and in the MO and HK regions during summer through increased density stratification and changes in water column mixing (Fig. 12). Increased river discharge is conducive to seaside transport in the UPPER and WEST regions during summer and in the HK region during winter, which is associated with an intensified offshore current and the seaward movement of the density front. The removal of tidal currents and reduced river discharge lead to an intensified landward current and westward transport current from Hong Kong waters and the adjacent shelf, accelerating bottom-water intrusion from the lower part to the upper part of the estuary.