The model used in this study was a fully validated hydrodynamic numerical Estuarine, Coastal and Ocean Model with a semi-implicit scheme (ECOM-si). It was modified from the Princeton Ocean Model, which was originally developed by Blumberg (1994). ECOM-si couples the sediment dynamics module and ecosystem module to simulate phytoplankton dynamics.

The sediment dynamics module combines suspended sediment transport equations with the hydrodynamic model. The Simulating Waves Nearshore wave module (Booij et al., 1999) was introduced to calculate the bottom shear stress induced by wave processes. Comprehensive sediment dynamics, including resuspension, deposition, settling were considered. The coupled sediment dynamics module has been successfully applied in the Changjiang River Estuary by Luo et al. (2017).

The main component of the ecosystem module used in this study was the vertical one-dimensional N2P2ZD model, which was developed from the Flexible Biological Module of the Finite-Volume Coastal Ocean Model (Chen et al., 2006) and was introduced into the EOCM-si. In recent years, the parameters of the ecosystem module have been further adjusted based on the variation and adjustment of phytoplankton in the areas off the Changjiang River Estuary, combined with in-situ observations in previous years (Wang et al., 2019a). In addition, the parameterization scheme for the effect of light availability on phytoplankton growth has been optimized, and has also been successfully applied on the phytoplankton dynamics in the Changjiang River Estuary (Wang et al., 2019a). The simulated components included nutrients (dissolved inorganic nitrogen and phosphate), phytoplankton (dinoflagellates and diatoms), zooplankton, and detritus as well as their photosynthesis, respiration, mortality, and grazing activities under specific conditions. The model performance of phytoplankton and nutrients has been well validated in our previous studies (Wang et al., 2019a; He et al., 2022; Wang et al., 2023), exhibiting reasonable spatiotemporal distribution of phytoplankton and nutrients.

The computational domain was roughly fan-shaped and covered the entire Bohai Sea, Yellow Sea, and East China Sea as well as part of the Japan Sea and the western Pacific Ocean ( Figure 1A ). The model grid had a high resolution in the horizontal direction, reaching about 1 km in the Changjiang River Estuary and 3–4 km near the coasts. A sigma coordinate of 20 layers in the vertical direction was used to fit the topography. The hydrodynamic open boundary was driven by the momentum flux, which includes tidal currents and shelf currents. The tidal current boundary was defined by 11 subtidal harmonic constants, while the shelf current boundary was derived from the monthly Simple Ocean Data Assimilation (SODA) datasets. The open boundary conditions and initial temperature and salinity values were also derived from the same datasets. The sea surface boundary heat flux was calculated based on sea surface temperature and atmospheric parameters (Ahsan and Blumberg, 1999), including mean sea level pressure, 2-m dewpoint temperature, cloud cover, 2-m air temperature, and 10-m wind, which were obtained from the ERA-5 database of the European Centre for Medium-Range Weather Forecasts. The Changjiang River discharge data were obtained from the Datong hydrological station, which is located 630 km upstream of the estuary. This hydrodynamic model has been successfully applied to simulate the dynamics of the Changjiang diluted water (Wu et al., 2011; Yuan et al., 2016; Wu and Wu, 2018; Zhang et al., 2018). The data reported by Yang et al. (2015) were used to set the upstream suspended sediment flux boundaries, with the initial and open boundary conditions of the suspended sediment being set to zero. The upstream nutrient flux boundaries were set based on Gao et al. (2012). And the nutrients’ initial and open boundary conditions were derived from the climatological monthly data of the World Ocean Atlas (WOA09). The atmospheric deposition of nutrient was not considered, since the flux in the coastal region off the Changjiang River Estuary was negligible compared to the eutrophic condition resulting abundant fluvial nutrient load. The initial conditions for other ecological variables in the model were set to be uniform background values and so were the open boundary conditions of these ecological variables. The model was spin-up for two years to achieve stable and satisfactory results.

The dynamic conditions of the Changjiang River Estuary and its adjacent waters are very complex and are subject to the combined effects of runoff, tides, wind, and other factors. Therefore, numerous related sensitivity experiments were designed to investigate the effects of key physical processes on the distribution of phytoplankton over multiple time scales. All conditions of the control experiment (Case 0) and other conditions of the sensitivity experiments were climatological. All sensitivity experiments with changed conditions were hot-run in the third year (full-year run started from January or flood-season run started from July). The details of the model settings are shown in Table 1 .

Case 1 was designed without tidal forcing, maintaining all the other conditions consistent with control experiment. The results of this case could give us a background variations of chlorophyll without tide.

Case 2 was designed with a 7-day tide delay compared to the control experiment, while maintaining all the other conditions consistent with control experiment, to distinguish distribution and variation such as phytoplankton and nutrients between cases during spring and neap tides (Case 0 and Case 2).

The Changjiang River carries a large amount of freshwater and terrestrial materials, which have a significant impact on the outbreak of phytoplankton blooms in the estuary area and its adjacent waters. Numerous extreme weather events have affected the Changjiang River in recent years, resulting in either extremely high or low discharges ( Figure 2A ), specifically in the cases of large floods and droughts, respectively. Such differences in the amounts of freshwater and nutrients have influenced the hydrodynamic conditions of the study area to different degrees, with varying effects on the intensity and scale of phytoplankton blooms. Therefore, two more sensitivity experiments, Case 3 and Case 4, were set up, in which the Changjiang River discharge increased and decreased by 50%, respectively, based on the control experiment.

Atmospheric conditions play an essential role in the occurrence of phytoplankton blooms. Wind speed and wind direction have not been stable in recent years, and the northerly winds in summer have not been accidental ( Figures 2B–E ). Southeasterly winds are predominant for most of July in the Changjiang River Estuary and its adjacent waters in summer, except for a small area where southwesterly winds prevail. As global warming increases, typhoons are occurring more frequently and are characterized by higher wind speeds (Bhatia et al., 2018). Changes in the wind field are often accompanied by variations in the vertical mixing of the water column, which modify nutrient distribution patterns and other conditions, thus indirectly affecting the growth and distribution of phytoplankton and the occurrence of phytoplankton blooms. Therefore, in the next two experiments, Case 5 and Case 6, different wind speeds during the flood season were simulated to investigate the effect of summer wind speed on phytoplankton blooms.

A last set of experiments, Case 7 to Case 20, along with Case 3 and Case 5 which mentioned previously were designed to establish the effects of the time scales of dynamic processes on phytoplankton blooms.

In this study, the “probability of phytoplankton blooms” was used as an index to describe the scale and degree of the chlorophyll maximum in a specified period. The probability of occurrence was calculated as follows:

where A _i represents the serial number of the statistical unit in the study area; T 1 and T 2 are the start and end times of the participation statistics, respectively; and p ( A i ) is the occurrence probability of phytoplankton blooms, which was defined by the following equation:

where C h l t h r e s h o l d is the defined high chlorophyll concentration threshold, i.e., 5 mg/m³ (Korb et al., 2004; Sòria-Perpinyà et al., 2021).

The EOF was used for the analysis of the daily chlorophyll values obtained from the model to understand the intra-annual spatiotemporal variability of chlorophyll distribution. The EOF, also known as principal component analysis, is a common method applied in mathematical statistics to analyze random variables (Wei, 2007). It is particularly useful for the analysis of spatial and temporal patterns in elemental studies and is commonly used in marine and atmospheric sciences (Thompson and Wallace, 1998; Storch, 2001).

The key indices of the horizontal chlorophyll maximum in different sensitivity experiments (i.e., area proportion, mean value, and offshore distance) were calculated to determine the response time of phytoplankton distributions to physical processes. The area proportion was calculated as follows:

where ∫ p ( t ) d A i represents the area with a chlorophyll concentration greater than 5 mg/m³, t represents time, and A represents the study area. The mean value is the average chlorophyll in areas with chlorophyll concentrations greater than 5 mg/m³. The offshore distance is the average distance between area with chlorophyll concentrations greater than 5 mg/m³and position (122°E, 31°N).

The results of the numerical model were used to further analyze phytoplankton variability. Seasonal variations in chlorophyll were analyzed based on the monthly average surface chlorophyll concentration and the probability of phytoplankton blooms in representative months in the control experiment (Equations 1, 2; Figure 4 ). The results showed that the monthly surface chlorophyll values were high in spring and summer (represented by April and July, respectively, Figures 4B, C ), decreased in autumn (represented by October, Figure 4D ), and remained low in winter (represented by January, Figure 4A ). Phytoplankton blooms started to appear sporadically in April, peaked in July, and then started to decrease, receding and shrinking in October, and almost disappeared in January. The horizontal chlorophyll maximum region was located between 122°E and 123°E, parallel to the coastline, and appeared as a long, curved, and oscillating strip with small patches from the northeastern to the southern end of the Changjiang River Estuary. The results obtained from the numerical model showed that the bloom distribution patterns were consistent with the phytoplankton bloom records.

It should be notable that our numerical model results showed the highest biomass in July ( Figure 4 ) while the results of phytoplankton bloom events ( Figure 3 ) showed the highest frequency in May and June. This discrepancy was largely due to the seasonal succession of predominate algae species in the study region (Tang et al., 2006; Yang et al., 2014). The complex biological processes including grazing and species competitions increased the simulation difficulty. Even so, previous field investigations demonstrated that the highest biomass emerged in summer season in the study region (Yang et al., 2014). The discrepancy was beyond the main purpose of this study.

Indeed, besides the seasonal variations, the horizontal chlorophyll maximum region also has variations within monthly scale. The monthly probability of phytoplankton blooms computed by the model showed that the horizontal chlorophyll maximum region varied significantly over periods shorter than one month ( Figure 4 ). For instance, some identified horizontal chlorophyll maximum region in April only had the value of 50% probability of phytoplankton blooms ( Figure 4B ). Even in July, the core region of probability could not reach 100% ( Figure 4C ), implying that the horizontal chlorophyll maximum region varies at the spatiotemporal scale, i.e., by oscillating in space, growing and declining over short time scales.

The daily averaged chlorophyll time series at five sites (shown in Figure 1B ) in control experiment (Case 0) were selected to further analyze the variation of phytoplankton biomasses in the Changjiang River Estuary at more specific time scales. These sites were arranged either across or along the shelf to represent spatial differences in both directions ( Figure 1B ). The hourly chlorophyll data at each site from April to September in Case 0 and Case 1 were also spectroscopically analyzed to determine the main time scale of phytoplankton variation.

In general, the spectroscopic analysis of all five sites in control experiment detected periods of approximate half a day, one week, and one fortnight, corresponding to the signals of flood-ebb, half-spring-neap, and spring-neap cycles (red lines in Figures 5F–J ). However, the results with tidal closure showed significantly weakened periodic signals (sky blue lines in Figures 5F–J ). Only periods of shorter than half a day, and one week was detected, while the period of one fortnight was missing. These periodic signals of chlorophyll variations without tide were corresponded to the signals of diel solar irradiation and inherent periodic signals of chlorophyll variations. The daily averaged chlorophyll time series showed that the periodic signals of approximately one fortnight owned large amplitudes, while those of one week presented small fluctuations ( Figures 5A–E , also seen in Supplementary Figure S1 ). This implied that the spring-neap tidal variability dominated phytoplankton variations over the time scale of days in spite of some inherent minor fluctuations of the half-spring-neap cycles.

The time-series data at sites 1 to 3 showed the differences in phytoplankton biomass across the shelf. These positions also crossed the horizontal chlorophyll maximum region east of the Changjiang River mouth, an area where both onshore and offshore variations in phytoplankton biomass can be observed. In spring and summer, the data at site 1 showed significant high value of phytoplankton with obvious variations (grey shade), while that at site 2 showed longer phytoplankton blooms, indicating that more distinct event-scale variations occurred onshore than offshore and more distinct phytoplankton blooms occurred offshore than onshore in the horizontal chlorophyll maximum regions ( Figures 5A, B ). That at site 3 showed the smallest variation and seldom high value of phytoplankton ( Figure 5C ), indicating that phytoplankton seldom blooms in this region. Moreover, the results of spectroscopic analysis demonstrated that the period signal was stronger onshore than offshore in the regions of horizontal chlorophyll maximum. The above results implied that 1) in the regions of the chlorophyll maximum, phytoplankton blooms statistically lasted longer offshore than onshore, and 2) phytoplankton blooms onshore occurred relatively periodically while those offshore occurred randomly ( Figures 5F, G ).

The spectroscopic analysis of site 2, 4, and 5 revealed the differences along the shelf. At these sites, from north to south, were also located small patches of the chlorophyll maximum region. The results showed that the daily variance signals became stronger when chlorophyll value sustained in high levels at the north two sites (site 2 and site 4). For instance, the shaded area was larger in June and July than in April and May at both Points 2 and 4; however, the signal was not significant at Point 5 ( Figures 5B, D, E ). The above results suggested that the chlorophyll maximum regions that were close to the Changjiang River Estuary exhibited more significant event-scale variabilities. Furthermore, variations occurred over short time scales in spring and summer, which was also consistent with the results shown in Figure 3C indicating that the durations of phytoplankton bloom events were last for days from May to July.

Both the records of the phytoplankton bloom events and the numerical model results demonstrated that phytoplankton biomass varied spatially over multiple time scales, from the seasonal to the weekly and daily scales. To reveal the regularity of chlorophyll variations over multiple time scales in the whole study area, the EOF was then used for further analysis of the numerical model results. The analysis was based on the daily averaged chlorophyll concentrations of the model, and three main EOF modes were calculated, each of them accounting for 47.3%, 17.7%, and 7.9% of the total chlorophyll variation rate, i.e., 72.9% of the total variation rate ( Figure 6 ). This meant that the three main modes could represent the chlorophyll variation in the Changjiang River Estuary and its adjacent waters. The area with high standard deviation of chlorophyll was mainly located between the 30-m and 50-m isobaths, which corresponded to the usual location of the horizontal chlorophyll maximum region.

The variance contribution rate of mode I was 47.3%, which represented the most dominant trend of chlorophyll variation and reflected the main characteristics of the seasonal variation of phytoplankton. Spatially, the horizontal chlorophyll maximum region extended to 30°N and south of the area. Based on the variation of the temporal coefficients, this region was characterized by significant seasonal variation. The horizontal chlorophyll maximum region outside the Changjiang River Estuary in spring and summer was the most dominant distribution characteristic of phytoplankton in this area. This region began to form in May and lasted until the end of September. However, the temporal coefficients also fluctuated over the time scale of weeks.

The variance contribution rate of mode II was 17.7%. This mode described the evolution of the horizontal chlorophyll maximum region, from generation to recession. The phytoplankton biomass first grew along the Zhejiang coast and the shelf area outside the Changjiang River Estuary in April and gradually gathered northward along the coast to form a horizontal chlorophyll maximum region, which started to decline toward the shelf and southward in August.

The variance contribution rate of mode III was 7.9%. This mode mainly showed the strong spatial oscillation of the horizontal chlorophyll maximum region over the time scale of days. It also showed that these regions extend onshore or offshore, and the intensity of these regions decreased or increased over this time scale. The temporal coefficients of modes II and III were also shown to be more likely to change its positive-negative sign in summer than in spring, indicating that phytoplankton was more likely to vary over short time scales in summer than in spring.

Both the time series of the modeled chlorophyll concentrations and the results of EOF analysis showed that the phytoplankton distributions presented high-frequency and stable oscillating signals ( Figures 5 , 6 ). The results of spectroscopic analysis also showed that the stable signals had periods of approximately half a day and a fortnight ( Figures 5F–J ), which corresponded to the flood-ebb and spring-neap tidal cycles in the study area. Thus, the frequency of these chlorophyll variations should be strongly influenced by tidal forcing.

The first two sensitivity experiments (Case 2 and Case 0) examined chlorophyll variability during the spring-neap tidal cycles. The results showed that the phytoplankton bloom probability was generally higher in Case 2 than in Case 0 ( Figure 7D ), specifically on the onshore and offshore sides of the chlorophyll maximum regions. These results indicated that phytoplankton blooms during spring tide probably occurred in a smaller region compared to those during neap tide. As the tidal cycle periods were relatively stable, the phytoplankton variations during different tidal cycles could be regarded as intrinsic, which were assisted in exploring the effect of other dynamic processes.

The runoff of the Changjiang River is a key dynamic process controlling the hydrological and environmental conditions in the estuarine region. The second and third sensitivity experiments (Case 3 and Case 4) examined chlorophyll variability under the increase and decrease of runoff during the flood season, respectively. The results showed that the increase and decrease of the Changjiang River discharge during flood season could produce opposite impacts on phytoplankton variation ( Figures 7B, C ). In the case of increasing runoff (Case 3), the probability of phytoplankton bloom occurrences largely increased on the offshore side of the chlorophyll maximum region and decreased on the onshore side. The opposite was observed in the case of decreasing runoff (Case 4).

A number of spatial differences in the probability of bloom events due to the effect of runoff were also detected, and they were more significant in the northern parts than in the southern parts of the study area. For instance, the variant values of probability were higher in the region around site 2 and 4, and lower in the region around site 5. These results suggested that the intensity of influenced area by runoff was gradually weakened as the distance prolonged from the river source. Furthermore, the influence on the southern parts adjacent to the Changjiang River Estuary was also weaker than that on the northern parts during summer season.

Wind is another key element affecting dynamic processes in riverine environments. The fourth and fifth sensitivity experiments (Case 5 and Case 6) examined chlorophyll variability under the increase and decrease of the prevailing summer monsoon, respectively. The results showed that under the condition of increased summer monsoon ( Figure 7E ), the probability of phytoplankton blooms increased in the core and offshore side of the chlorophyll maximum region, while it decreased on the onshore side. However, the decrease in value was smaller than the increase in value. This implied that the enhanced wind conditions stimulated the occurrence of phytoplankton blooms and favored the shifting of the chlorophyll maximum region offshore. In contrast, as the summer monsoon declined ( Figure 7F ), the probability of phytoplankton blooms at the core of the chlorophyll maximum region decreased. However, the affected region in this case was not as large as that under the enhanced wind conditions ( Figure 7E ), specifically at the southern edge of the study area. The region with an increased probability of phytoplankton blooms on the onshore side of the chlorophyll maximum region in Case 6 was also not obvious, which implied that the reduced wind conditions were not conducive to the maintenance of phytoplankton blooms.

When the runoff pulses increased during spring tide, the response of area proportion took more than 10 days ( Figure 8A ); in contrast, when they increased during neap tide, the response occurred in less than 5 days ( Figure 8B ). These sensitivity experiments (Case 3, Case 7-Case 13) involving increased runoff pulses showed that the responses of area proportion were always negative during neap tide when the horizontal chlorophyll maximum developed. These responses gradually shifted and became positive when the horizontal chlorophyll maximum receded. This situation only happened when the runoff pulses could last until the phytoplankton blooms developed from spring to neap tides. Otherwise, the negative response would gradually disappear or even never emerge.

The response of mean value was different from that of area proportion ( Figures 8C, D ). In all the sensitivity experiments with increased runoff pulses, the mean values showed a slightly negative response during the first 3 days. After this time, the mean values began to show a positive response, which lasted up to 2–4 days after the runoff pulses ended. Specifically, this delayed positive response was slightly stronger than that occurring during the runoff pulses.

The response of offshore distance to runoff was always positive. All the sensitivity experiments showed an integrated increased response of offshore distance ( Figures 8E, F ). The variations of offshore distance gradually recovered when the runoff pluses ended. Indeed, the offshore distance would not persistently increase, it maintained a relatively stable value after 7–15 days of runoff pluses even if the pluses still sustained.

When wind pulses increased during spring tide 5 days were required for the response of area proportion ( Figure 9A ). In contrast, when they started during neap tide, the response took less than 5 days ( Figure 9B ). The area proportion showed mostly positive responses during neap tide when it increased and reached its maximum. This occurred even when the wind pulses ended before the maximum area proportion was reached; for instance, in the case of wind pulses lasting for over 7 days during spring tide and over 3 days during neap tide.

The response of mean value was more rapid than that of area proportion ( Figures 9C, D ). The results showed immediate positive responses of mean value within the first 2 days. Subsequently, it slowed down and fluctuated, and remained stable for about 2 weeks if the wind pulses persisted. Once these ended, the positive response of the mean value would also recover to normal in 1 or 2 days.

The response of offshore distance was also mono-positive. The rapidly increasing response of offshore distance occurred in 4–6 days and could be maintained for 5 more days if the wind pulses persisted ( Figures 9E, F ), but it decreased when the wind pluses ended. The periods of rapidly increasing response also matched the timing of the horizontal chlorophyll maximum extending offshore. Offshore distance also did not persistently increase; it remained relatively stable during wind pulses (for 15–22 days) and even after they ended.

Although the increased runoff and wind caused similar variation patterns in the horizontal chlorophyll distribution, differences were observed in relation to the time scale. For example, the effect of wind was faster than that of runoff ( Figures 8 and 9 ). All of the three indices (area proportion, mean value, and offshore distance) of the horizontal chlorophyll maximum responded more rapidly to wind than to runoff. The responses to wind pulses occurred a few days earlier and were stronger than that to runoff. At the same time, the periods of recovery were shorter when the wind pulses ended than when the runoff pluses ended. Another difference was that delays were observed more for the runoff effects than for the wind effects. For instance, the response of mean value was noticeably delayed when the runoff pulses were over ( Figures 8C, D ), but such delay was not observed when the wind pulses ended ( Figures 9C, D ).