Long-term response of an estuarine ecosystem to drastic nutrients changes in the Changjiang River during the last 59 years: A modeling perspective

Shi, Shenyang; Xu, Yi; Li, Weiqi; Ge, Jianzhong

doi:10.3389/fmars.2022.1012127

ORIGINAL RESEARCH article

Front. Mar. Sci., 17 October 2022

Sec. Marine Biogeochemistry

Volume 9 - 2022 | https://doi.org/10.3389/fmars.2022.1012127

This article is part of the Research Topic Multi-scale Variability of Ecosystem Functioning in European and Chinese Shelf Seas View all 16 articles

Long-term response of an estuarine ecosystem to drastic nutrients changes in the Changjiang River during the last 59 years: A modeling perspective

Shenyang Shi¹

Yi Xu^2*

Weiqi Li²

Jianzhong Ge^2,3

¹Project Management Office of China National Scientific Seafloor Observatory, Tongji University, Shanghai, China
²State Key Laboratory of Estuarine and Coastal Research, East China Normal University, Shanghai, China
³Institute of Eco-Chongming, Shanghai, China

The riverine nutrient inputs to the ocean reflects land-use changes and can affect the health of coastal environments over time, especially for a highly-anthropogenically influenced river-estuary-shelf system. To investigate the impact of riverine inputs on the Changjiang Estuary ecosystem at a multi-decadal time scale where long-term observations are limited, we built a three-dimensional physics-biogeochemistry-coupled model system based on the Finite-Volume Community Ocean Model (FVCOM) and the European Regional Shelf Ecosystem Model (ERSEM). Our model successfully simulated the temporal and spatial nutrient variabilities in the river-estuary-shelf con7tinuum from 1960 to 2018. The results showed increasing trends of nitrate and phosphate and fluctuating silicate variability, thereby leading to rising nitrogen (N) to phosphorus (P) ratios and decreasing silicon (Si) to N and P ratios. Such changes in the stoichiometric relationship of nutrient species also alter the community structure of the primary producers in estuaries. Our model showed a general increase of diatoms over the 59 years, corresponding to decreased proportions of micro-phytoplankton and pico- phytoplankton. With different backgrounds of light and nutrient limitations in the river and inner shelf, our model suggests that the trend of the diatom proportion in the light-limited river mouth is more associated with silicate variability, with decreased diatom proportions occurring in the 2000s. Our model relates the hydroclimate, nutrient load, and biogeochemical cycling, reproducing estuarine ecosystem variability and clarifying issues such as the causality of the ecosystem interactions.

Introduction

Human civilization originated in areas around large rivers; the ecological environment of rivers and estuaries is closely related to human life. Agriculture and industry developed rapidly after the Second World War, leading to extensive use of fertilizers and detergents (Bennett et al., 2001; Galloway and Cowling, 2002; Duan et al., 2007). Disposal of sewage and wastewater has caused the water quality of rivers and estuaries to deteriorate. While pollution accumulates in rivers, estuaries are influenced by both point and non-point pollution sources, which has changed such ecosystems significantly in the last century (Rabalais et al., 2009; Wang et al., 2013; Powers et al., 2015; Dai et al., 2016). Large rivers are essential in this terrestrial-estuary linkage. Therefore, anthropogenic perturbations inevitably affect rivers, and so are estuaries (Howarth et al., 2000; Scavia et al., 2000; Howarth and Marino, 2006; Bricker et al., 2008; Oviatt et al., 2017). To some degree, nutrients and phytoplankton can be indicators of productivity and ecological health of estuaries (Glé et al., 2008; Enawgaw and Lemma, 2018). Over the past few decades, many large estuaries around the world have experienced dramatic ecological changes. In North America, the Mississippi River has shown a drastic increase in nitrogen and phosphorus since the 1950s, with dissolved inorganic nitrogen (DIN) and total phosphorus frequently exceeding 150 μM and 7 μM, respectively during the 1980s (Lohrenz et al., 2008), and riverine silicon has decreased since 1960 because of the stimulation and burial of freshwater diatoms (Rabalais et al., 1996; Houser et al., 2010; Leong et al., 2014). Therefore, the nutrient ratio and phytoplankton species composition in the coastal area of the Mississippi River has changed (Dortch et al., 2001; Smits et al., 2019; Stackpoole et al., 2021). In Europe, the nutrient concentration in the Scheldt Estuary increased significantly from 1965 to the mid-1970s and decreased from 1980 to 2002, which was positively correlated with the river discharge. Because of the difference in nutrient change ratios, the function of the estuary has shifted from a net nitrate producer to a net consumer (Soetaert et al., 2006; Middelburg et al., 2011).

The Changjiang Basin is 1.8 million square kilometers in total, accounting for 18.8% of China’s land area (Zheng et al., 2020). It has abundant natural resources and is the most densely populated area in China. Like other giant rivers and estuaries, anthropogenic activities significantly affect the Changjiang Basin and estuary (Huang et al., 2001; Chen et al., 2015; Xu et al., 2015; Dai et al., 2016). Many studies have focused on ecosystem changes in the Changjiang Estuary, including eutrophication, nutrient stoichiometry, and a regime shift. Hydraulic engineering and wastewater discharge are considered as a critical factor impacting the Changjiang Estuary ecosystem during the past few years (Yang et al., 2006; Gao and Wang, 2008; Dai et al., 2011; Li et al., 2015). As of 2018, over 50,000 dams have been built in the Changjiang Basin. These dams mediate water volume and mitigate floods while decreasing sediment discharge to approximately 140 Mt/y (Dai et al., 2008; Hu et al., 2009; Yang et al., 2018b). Concurrently, dissolved silicate has decreased by 30% from the 1960s to the 2000s due to large amounts of sediment trapping by the dams (Dai et al., 2011; Wang et al., 2018). Agriculture and industry began to develop in the 1960s, drastically increasing nitrate and phosphate loads which are associated with fertilizer usage and domestic sewage (Xing and Zhu, 2002; Duan et al., 2007; Zhou et al., 2008; Yan et al., 2010; Liang and Xian, 2018). From the 1960s to the 1990s, the amount of nitrogen fertilizer rapidly increased 15-fold from 1.0×10⁸ kg to 1.50×10⁹ kg, and phosphorus fertilizer increased approximately 12-fold from 4×10⁷ kg to 5.0×10⁹ kg (Duan et al., 2000; Duan et al., 2007; Shen and Liu, 2009; Shen et al., 2020). The use of nitrogen and phosphorus fertilizer still increased steadily after 1990s (Liu et al., 2022). As a result, the DIN concentration increased threefold, and the concentrations of dissolved inorganic phosphorus (DIP) and dissolved organic phosphorus increased by 30% (Yu et al., 2012; Liu, 2017). Because of the different nutrient loading rates, the N: P: Si ratio has changed (Jiang et al., 2014). For example, the N/P ratio was ~30-40 in the early 1960s, and it was reported to reach 103 and 268 in 1985 and 1997, respectively (Gu et al., 1981; Wong et al., 1998; Shen and Liu, 2009); a factor of 17-fold greater than the Redfield ratio (Si:N:P: ~16:16:1). The increased input of N and P into an estuary can result in eutrophication, directly or indirectly linked to harmful algal blooms (Wang, 2006; Shen and Liu, 2009) and hypoxia (Wang et al., 2016a). Concerning the bloom frequency in the Changjiang Estuary, diatom blooms have decreased due to silicate trapping. Dinoflagellate blooms tend to be more frequent because of nutrient increases (Jiang et al., 2014).

Our knowledge of the Changjiang nutrient flux is essentially based on hydrographic data measured at river gauge stations, which are usually located upstream (Zhang et al., 2021). How these riverine nutrients are processed downstream and their contribution to the estuary ecosystem over a longer time scale are poorly understood. The times and places of the existing studies of Changjiang Estuary are limited (Li et al., 2007; Zhou et al., 2008; Jiang et al., 2014; Yang et al., 2015; Luan et al., 2016; Jiang et al., 2017; Jiang et al., 2018a). Available studies mostly used sediment core sampling and satellite datasets (Henson et al., 2010; Jin et al., 2010; Cheng et al., 2014; Zhu et al., 2014; Chen et al., 2017; Zhang et al., 2021). Therefore, the data based on riverine hydrographic stations and dispersed sampling could be biased concerning the seaward river fluxes and the ratios of nutrient species related to estuary ecosystem dynamics (Zhang et al., 2021).

Numerical models coupled with physical and biogeochemical dynamics are an effective way to improve our understanding of the impact of riverine nutrient inputs on ecosystem function in this land-sea boundary. Different models have been developed to compensate for the lack of empirical data and study how ecosystems respond to environmental changes (Riley, 1949; Fasham et al., 1990; Sarmiento et al., 1993; Franks and Chen, 2001; Spitz, 2003; Gruber et al., 2006). In Ge et al. (2020a), the Finite-Volume Community Ocean Model (FVCOM) and the European Regional Shelf Ecosystem Model (ERSEM), was used to study the effect of offshore sediment front on seasonal marine ecosystem dynamics; the interannual variability of nutrient and phytoplankton in the Changjiang Estuary and shelf from 1999 to 2017 were discussed based on this model in Ge et al. (2020a). However, none of these studies involved simulations on a multi-decadal scale. In this study, based on a delicate grid design that can refine different spatial scale from the river channel to the estuary and the shelf, we used the same FVCOM and ERSEM coupled model to simulate long-term (59 years) variations of major ecosystem components along the Changjiang Estuary, with the objectives of (1) presenting the historical state, including the temporal and spatial distribution of nutrients in the river-estuary-shelf continuum; (2) investigating when and where the phytoplankton communities have changed in the estuary and inner shelf and how these changes relate to the trends of nutrient variability. For this study, we hope to use the coupled model driven by long-term empirical upstream nutrient data to assess the present state of the Changjiang Estuary and compare it to previous states.

Data and methods

Study domain and long-term evolution of river inputs

Our study area included the downstream of the Changjiang Basin, west of Datong (DT) station, to the east at 124°E, representing a continuum of river channel, river mouth, and inner shelf (Figure 1). Regarding river discharge, we only considered the Changjiang River, the fourth largest river in the world (Milliman and Farnsworth, 2011), originates in the Qinghai-Tibet Plateau (Wang et al., 2016b), flows through 19 provinces, and then into the East China Sea (Figure 1A). Until the end of the twentieth century, DT (Figure 1A) was still the most downstream hydrographic station in the Changjiang Basin. River discharge used in this study was collected at DT (Figure 1A) (http://xy.cjh.com.cn/). Nutrient data at the DT station and nearby sites were collected from different sources, i.e., previous studies conducted from 1960 to 2019 (Dai et al., 2011; Gao et al., 2012; Zhu et al., 2014; Li et al., 2021), and our own measurements made since 2013 (Ge et al., 2020a; Ge et al., 2020b; Ge et al., 2021). A detailed description of the nutrient sources can be found in Table 1.

FIGURE 1

Figure 1 Diagram of the study area [modified from Ge et al. (2020a)] (A) Location of the Changjiang Basin, the Three Gorges Dam (TGD), Datong Station, and Xuliujing Station. (B) Changjiang Estuary and inner shelf of East China Sea. Site 1 (at the river channel), site 2 (at the river mouth), and site 3 (at the inner shelf) are typical sites that were chosen (red dots). The dashed line denotes the boundary of the turbidity maximum zone.

TABLE 1

Table 1 Different sources of nutrient data collected from 1960 to 2018.

We collected approximately 60 years of annual mean surface runoff and riverine nutrient data (Figure 2). Significant interannual variability occurs along this time series. The strongest inundation in 1998 corresponded with the most extensive runoff in 1998 in the time series (over 1.2×10¹² m³/yr). In contrast to the fluctuating annual discharge, the sediment load clearly decreased, from to 500 Mt/yr to 100 Mt/yr continuously (the red curve in Figure 2A). This phenomenon is consistent with the results in Luan et al. (2016), that sediment was trapped in reservoirs mainly due to the construction of water conservation projects (Yang et al., 2015; Luan et al., 2016; Yang et al., 2018b). Riverine nutrient variabilities also showed distinct trends. Increasing trends from 1960 to 1980 were clearly seen in the nitrate and phosphate time series (with slopes of 1.98 (p<0.001) and 0.02 (p<0.001), respectively). By 2018, the nitrate and phosphate levels had increased by approximately 15-fold and two-fold, respectively, compared to the early 1960s. In contrast, silicate showed only a slightly increasing trend (Figure 2D). The Pearson correlation coefficients for these three nutrients were less than 0.001. Therefore, the fitted curves showed significant correlation, and the trends for the nitrate, phosphate, and silicate concentrations were reliable. However, there is little empirical data for ammonium, so that trend line does not accurately represent the ammonium concentration (Figure 2E).

FIGURE 2

Figure 2 Empirical annual water runoff data (blue line) and sediment load (red line) at Datong Station (A). Observed (red dots) and empirical data (blue lines) for nitrate (B), phosphate (C), silicate (D), and ammonium (E) at Datong and nearby sites.

Coupled model system and configuration

To describe the long-term variations of ecosystems in the Changjiang Estuary continuum during dramatic environmental changes, we used a physics-biogeochemistry coupled FVCOM-ERSEM model to describe the variability of major nutrients and the phytoplankton community during the past few decades. FVCOM (Chen et al., 2003) discretizes space into an unstructured triangular grid in the horizontal direction and uses terrain tracking coordinates vertically; thus, it is suitable for the complex situations that characterize estuarine and coastal areas. The use of scalar conservation equations makes FVCOM user-friendly for multidisciplinary research. The ERSEM model is one of the most comprehensive models for the low-trophic marine food web (Butenschon et al., 2016) and has been widely used in the simulation of biogeochemical processes in coastal areas (Elkalay et al., 2012; Sankar et al., 2018; Jardine et al., 2021). The ERSEM has complex ecosystem components, including nutrients, phytoplankton, zooplankton, benthos, bacteria, iron, light extinction, calcification, and alkalinity. In ERSEM, four functional types of primary producers were considered (diatoms, micro-phytoplankton, pico-phytoplankton, and nano-phytoplankton), which allowed us to identify the long-term changes in phytoplankton communities when nutrient concentrations and ratios change. The coupling of these two models has been successfully applied in physical and ecosystem modeling in the Changjiang Estuary and adjacent regions for both episodic events and decade-scale modeling (Ge et al., 2020a; Ge et al., 2020b; Ge et al., 2021).

The unstructured triangle grid used in this study is of medium-resolution and covers the entire study area. It covers the Changjiang Estuary, Hangzhou Bay, and the inner shelf of the East China Sea. In the coastal region, river mouth, and in other complex terrains, the resolution of the grid was adjusted to 2-3 km to ensure that detailed hydrodynamic processes could be geometrically fitted in the model. The hydrodynamic model FVCOM involves FVCOM-SED and FVCOM-SWAVE with current-wave–sediment interactions (Wu et al., 2011). The simulation period was 1960–2018, and the output was daily averaged. The nutrients and dissolved oxygen of the open boundary were monthly climatologies obtained from the 2013 World Ocean Atlas dataset. The model was driven by a reanalysis of atmospheric forcing data at 3-h intervals, provided by the National Centers for Environmental Prediction from 1960 to 1980 and datasets of ERA-Interim from the European Center for Medium-Range Weather Forecasts (ECMWF) from 1980 to 2018. We used surface air temperature, pressure, relative humidity, winds at a 10m altitude above the sea surface, downward longwave radiation, and net shortwave radiation to specify the surface fluxes of momentum and buoyancy based on bulk formulae (Fairall et al., 2003). At the open boundary, we specified the temperature, salinity, nitrate, total inorganic carbon, alkalinity, and dissolved oxygen using the output of a regional model in this study domain described in Ge et al. (2013). The major tidal components along the lateral open boundary, including eight primary harmonic constituents M2, S2, N2, K2, K1, P1, O1, Q1, were determined from the TPXO 8-atlas (Egbert and Erofeeva, 2002). The riverine inputs of nutrients and freshwater were based on the empirical data plotted in Figure 2. In our simulation, we also included the benthic components of the ERSEM to consider the effects of the benthic ecosystem on the pelagic ecosystem. The flux from the benthos to the water column was simplified as a remineralization process from dissolved organic matter to inorganic matter.

The Empirical Orthogonal Function (EOF) analysis

EOF is a method that extracts an eigenvalue from a large-scale dataset. Such EOF analyses were widely used to compress the spatial and temporal variability of large datasets down to a set of dominant, mathematically orthogonal (independent) spatial functions and associated time-varying amplitudes (PCs) (Briggs, 2007; Dawson, 2016). It has been widely used in meteorology and oceanography (Lian and Chen, 2012; Qi et al., 2014; Roundy, 2015; Zhong et al., 2020; Grams, 2021). The EOF analysis was used to evaluate when and where the most dominant phytoplankton structure variabilities occurred during the 59-year time series.

Results

Verification of simulated nutrients and phytoplankton

To investigate nutrient and phytoplankton, we selected three specific sites downstream of Xuliujing to 124°E, which represented the river channel, the river mouth, and the inner shelf (red dots in Figure 1B) dynamics during the past 59 years. We used various in-situ observational data to verify the reliability of our numerical biogeochemical model results. We used various in-situ observational data to verify the reliability of our numerical biogeochemical model results. Based on previous studies, the nitrate concentration in the 1960s in Changjiang was ~20 μM and was maintained at this level for the next two decades. In the 1980s, the DIN concentration increased to 80 μM, then exceeded 100 μM in the 2000s (Jiang et al., 2010; Jiang et al., 2014; Liang and Xian, 2018; Shen et al., 2020). This is consistent with our simulated nitrate concentration at the site near Xuliujing station (Figure 3A, site 1). For the saline sites, the documented nitrate concentration was ~10 μM in the 1960s, which increased to ~30 μM in the 1980s, then continued to grow to ~45 μM in the 2000s (Jiang et al., 2014). Our simulated time series of nitrate concentrations at representative stations in the saline section (Figure 3B) exhibited the same magnitude as Jiang’s, however, it was slightly lower than that reported by Liang and Xian (2018). The observed phosphate concentration in the river site fluctuated from ~0.5 μM before 1980, increasing to ~1.5 μM in the 1990s (Shen et al., 2020). Our simulated phosphate (Figure 4A) are larger than those in Jiang et al. (2010), which were ~0.5 μM before the 1980s, ~0.7 μM during the 1980s and ~1.2 μM in the 2000s. After 2003, the phosphate concentration increased quickly with sharp fluctuations (Figure 4A). The fluctuations since 2003 can be found in Liang and Xian (2018); however, they occurred early in our simulation. For the phosphate concentration in the saline section (site 3), our model successfully reproduced the concentrations observed by Jiang et al. (2014) and Liang and Xian (2018) (Figure 4C). We didn’t see a significant trend for the simulated silicate concentration. It essentially remained between 100-150 μM from 1985 to 2014. This is consistent with the results of Jiang et al. (2014) and Liang and Xian (2018) (Figure 5A). For the saline site (site 3), the long-term silicate trend is similar to that of Jiang et al. (2014), however, the annual average silicate was less than 20 μM, but in Jiang et al. (2014), it fluctuated around 30 μM (Figure 5C). The simulated nutrient concentrations were in the same order of magnitude and indicated similar temporal variability as the observations.

FIGURE 3

Figure 3 Modeled time series of nitrate concentration at the three sites. (A: site1, B: site2, C: site3). Black lines demonstrate the daily averaged nitrate concentration, red lines are annual averaged values, and blue lines indicate ten-year moving averaged trend.

FIGURE 4

Figure 4 Modeled time series of phosphate concentration at the three sites. (A: site1, B: site2, C: site3). Black lines demonstrate the daily averaged phosphate concentration, red lines are annual averaged values, and blue lines indicate ten-year moving averaged trend.

FIGURE 5

Figure 5 Modeled time series of silicate concentration at the three sites. (A: site1, B: site2, C: site3). Black lines demonstrate the daily averaged silicate concentration, red lines are annual averaged values, and blue lines indicate ten-year moving averaged trend.

Data collection for the validation of the dominant phytoplankton species is more complicated. According to the sampling method, the empirical phytoplankton data can be divided into net-collected and water-collected phytoplankton. These different sampling methods can yield phytoplankton of different sizes and species. The parameters representing the condition of the phytoplankton, such as composition, abundance, and biomass, are also different. Phytoplankton also agglomerate; therefore, observations vary greatly over time and space. Hence, there was a deviation between the data from previous studies and the model results. Previous studies have documented that the annual mean diatom ratio in the phytoplankton communities was over 80% before 1990 and then decreased to 69.8% (Zhou et al., 2008; Jiang et al., 2014). The simulated data showed a similar trend, but the annual mean diatom ratio was lower, which varied in the range of 51.5% – 52% and 48.3 – 51.2% after 2000 at site2 and site3, respectively (Figure 7).

Simulated long-term variability of nutrients

Our model fit the observed data well for a specific time frame, it successfully simulated the seasonal and long-term nutrients variations at different geographical locations. Compared with site 1, which is influenced by riverine nutrients, the time series of nutrients at sites 2 and 3 are influenced more by oceanic water, thus indicating significant seasonal variations because of nutrient uptake by phytoplankton (Figures 3–5).

The nitrate concentration at the three sites showed a pronounced increasing trend from 1960 to 2019 (Figure 3). The Mann-Kendall test (Table 2) suggested a significant increasing trend (the τ value of the three sites was greater than 0.7, and the slope was positive with p< 0.001). Site 1, which is located in the river channel, had a nitrate concentration above 160 μM during the last decade. It was 140 μM for site 2, located in the river mouth, and much lower for site 3, located on the shelf. This suggests that the farther offshore the location, the lower the nutrient concentration. This phenomenon is associated with the geographical features in this continuum. The nutrient concentration is more greatly affected by freshwater dilution at the river mouth, while it is more influenced by the phytoplankton uptake offshore beyond the maximum turbidity zone. The same pattern can be observed from the ten-year moving average (blue lines in Figure 3).

TABLE 2

Table 2 Mann-Kendall test of nutrients and phytoplankton at three sites.

The long-term variability of nitrate during the 59 years of simulations can be divided into three phases: 1960–1980 (a steady increase), 1980–2010 (a rapid increase), and 2010–2018 (violent fluctuation). From 1960 to 1980, nitrate at the three sites was maintained at a relatively low level (~20 to ~30 μM in the estuary area and ~7 μM in the coastal region). The concentration gradually increased to 60 μM by the end of the 1970s. During the second phase, nitrate increased rapidly following the rapid economic growth in China since the early 1980s. The nitrate concentration at site 1 almost tripled owing to increasing emissions of industrial and agricultural sewage, which led to eutrophication in the Changjiang River during this period (Wang, 2006). During the last phase, the nitrate concentration fluctuated; however, it generally showed a steady increasing trend.

Similar to nitrate, the phosphate concentration increased at the three sites during the 59 years (Figure 4), but at a lower rate. Both site 1 and site 2 showed decreasing trends from 1960 to 1970, with the phosphate concentration falling from 0.7 μM to 0.1~0.2 μM. From 1970 to 2018, the concentration showed a steadily increasing pattern in the river channel (from ~0.2 μM to ~2 μM) and from ~0.2 μM to ~1.5 μM in the river mouth (Figures 4A, B). For the inner shelf at site 3, the phosphate concentration showed a slightly increasing trend compared with sites 1 and 2. From the Mann-Kendall test for phosphate, the slopes of these three sites were all positive (P< 0.001) (Table 2). The τ values for sites 1 and 2 were 0.67, but for site 3 it was 0.34, suggesting a weaker increasing trend at site 3 compared with the other two sites.

The temporal trends in the silicate concentrations suggested that silicate was the most stable nutrient component in the river-estuary-shelf continuum. There was no pronounced trend in the variability of silicate during our simulation period, especially for the river site (site 1) (Figure 5). This was also confirmed by the Mann-Kendall test, which suggested non-significant trends showing much higher P values (Table 2). Such a difference is reasonable because N, P, and Si are associated with different biogeochemical processes. The N and P concentrations have been strongly influenced by anthropogenic nutrient sources, whereas the primary silicate source within the Changjiang Basin is the weathering of rocks. From the time series for the silicate concentration, we can interpret the interannual fluctuations of silicate during specific periods. After 2003, the silicate content decreased significantly when the Three Gorges Dam (TGD) was completely put into operation (during 2003–2011: r = -0.46, p< 0.001), which suggests a retention effect of the TGD (Ding et al., 2019). For the oceanic site, a decreasing trend of silicate in the 2000s could be associated with its increased utilization by a diatom bloom as a result of increases in DIN and DIP after 2000.

Nutrient ratios and phytoplankton community composition

The nutrient stoichiometry also dramatically changed because of the different rates of change in nutrients. Generally, sites 2 and 3 showed similar increasing trends in the nitrate-to-phosphate ratio (N/P), with slopes of 0.38/a and 0.37/a, respectively (Figure 6A). The N/P ratio at site 2 was almost three times the Redfield ratio (16:1) by the end of 2018, and that of site 3 was close to the Redfield ratio before 1980 and increased significantly thereafter, reaching ~35:1 by the end of 2018. Our results suggest that water in the river mouth was more likely to be phosphate-limited during the 59 years, whereas on the inner shelf, there was no phosphorus limitation before 1980, but it became phosphate-limited after that.

FIGURE 6

Figure 6 Nutrient ratios of N/P (A), Si/N (B), Si/P (C) at site 2 (red) and site 3 (green). The dashed lines are the regression fitting of the nutrient ratios. Blue lines indicate threshold of the Redfield ratio.

Si/N and Si/P decreased owing to the rapid increase in N and P, while silicate evidenced a stable change (slope=-0.05/a and -0.03/a, respectively) (Figures 6B, C). At site 2, the Si/N ratio decreased from ~3.5-fold of the Redfield ratio (Si/N=1) to half of that at the end of our simulation. The cut-off point when the Si/N became less than 1 occurred around 2007, which suggests a possibility of silicate limitation since that time. The same variability occurred at site 3, but with a lower decreasing rate of Si/N. We attribute this to a less significant increasing trend of N at Site 3 when compared with the rapidly increasing N trend at site 2 (Figure 6B). The trend for Si/P also decreased at site 2; however, we did not find a significant decreasing trend at site 3 (Figure 6C), probably because of a much more drastic increase in the phosphate concentration at site 2 than at site 3.

Changes in the stoichiometric relationship of the nutrient species could alter the community structure of the primary producers in the Changjiang Estuary (Zhou et al., 2008; Chu et al., 2014). Previous studies have demonstrated that diatoms are the dominant species in the Changjiang Estuary (Jiang et al., 2014). Therefore, we assessed the variability of the diatom proportion in the phytoplankton community first. Our results were arrived at by calculating the percentage of diatom chlorophyll in the primary producers. In our simulation, diatoms dominate the phytoplankton community at sites 2 and 3, accounting for an average of 52.6% and 50.8%, respectively. As shown in Figure 7, the diatom proportions in the river mouth and inner shelf have been changed over the 59 years. Specifically, there was an increased frequency of diatom proportion by the year of the 2000, and decrease of that after. The long-term change of diatom proportion could be associated with changes in the riverine nutrient loads, including a consistent increase in nitrate and a decline in silicate after 2003 due to sediment trapping in the upstream reservoirs. With the decrease in Si/N after 2003, the proportion of diatom species began to decrease.

FIGURE 7

Figure 7 Time series of the diatoms proportion at (A) site 2 and (B) site 3. The black lines are ten-year filtered.

To further address the variabilities of other phytoplankton species on the estuary and inner-shelf. We checked the timeseries of micro-phytoplankton, nano-phytoplankton and pico-phytoplankton at site 2 and site 3 (Figure 8). It demonstrated that the fraction of micro-phytoplankton at site 2 (Figure 8A) is more than 20% and is larger than that of nano-phytoplankton (~15%) and pico-phytoplankton (~10%). However, at site 3 which located further offshore (Figure 8B), the nano-phytoplankton becomes the dominant specie among the three species. The proportion of micro-phytoplankton and pico-phytoplankton both remain at ~15%. As for the trended changes, the fraction of micro-phytoplankton at site 2 has declined since the late 1960s (Figure 8A). For site 3, both the pico-phytoplankton indicates the most significant decreasing trend (Figure 8B). These declines of smaller phytoplankton species correspond with the increase of diatoms proportions simultaneously, which is associated with the change of nutrient ratios.

FIGURE 8

Figure 8 Time series of other phytoplankton proportions at (A) site 2 and (B) site 3. P2 denotes nano-phytoplankton, P3 denotes pico-phytoplankton, P4 denotes micro-phytoplankton.

Spatial patterns of nutrients and phytoplankton variabilities

After analyzing the time series of nutrients and phytoplankton at the selected sites, we investigated how the changes in riverine nutrient inputs led to the temporal and spatial variability of nutrients and the associated phytoplankton distributions in this river-estuary-shelf continuum. We used the concentration in 1960 as the background value and calculated the percentage changes of nutrients in the next 10, 20, 30, 40, 50, and 60 years. From 1960 to 1969 (Figure 9), the nitrate concentration barely changed in the entire study area. As the nitrate input from the Changjiang increased gradually, it accumulated in the estuary and began to spread out to the sea. By the end of the 1970s, a small part of the inner shelf, the nitrate concentration increased by about 50% of its original value and nearly 80% in the river channel (Figure 9B). After 1980, the Changjiang Basin experienced a rapid nutrient increase over the next 30 years. In 1989, the river mouth and Hangzhou Bay nitrate concentrations were ~1.5-fold greater than those in 1960. Then, the nutrients started to spread further south of Hangzhou Bay (Figure 9C). This pattern persisted during the 1990s (Figure 9D). Until 2009 (Figure 9E), a large portion of the inner shelf was affected by nitrate from the Changjiang runoff, which could reach as far as 124.5°E in the east. At the river mouth, the nitrate concentration increased to 400% of the background value in 1960. During the 2010s, the nutrients from the Changjiang runoff remained stable because of wastewater discharge control (Shen et al., 2012; Shen et al., 2020). The riverine input affected area expanded north of the Jiangsu coastal area (Figure 9F). Phosphate showed a different pattern than nitrate (Figure 10), beginning to increase rapidly in the 2000s to a maximum of 200% (Figure 10E). At the same time, high phosphate levels also extended offshore, but unlike nitrate, this extension was limited to 123°E (Figures 10E, F).

FIGURE 9

Figure 9 Distribution of the changing percentage in the nitrate concentration during (A) 10 years, (B) 20 years, (C) 30 years, (D) 40 years, (E) 50 years, and (F) 60 years.

FIGURE 10

Figure 10 Distribution of changing percentage in the phosphate concentration during (A) 10 years, (B) 20 years, (C) 30 years, (D) 40 years, (E) 50 years, and (F) 60 years.

With variations in nutrient concentrations and nutrient proportions, the phytoplankton community also shifted in the whole estuary. While the nitrogen loads were accelerating, silicate did not increase significantly and even decreased after 2000; consequently, the Si/N ratios continued to decrease and have declined to lower than the Redfield ratio in recent decades. The silicate concentration is more likely to regulate diatom growth than other types of phytoplankton, e.g., dinoflagellates; thus, the diatom biomass decreased with the declining Si/N ratio. Therefore, the structure of the phytoplankton community changed. Our analyses of the diatom proportion over the whole study domain suggest that diatoms typically dominate the phytoplankton community of the Changjiang Estuary and inner shelf, contributing to more than 50% of the surface chlorophyll values (Figure 11A). To compress the spatial and temporal variability of the diatom proportion maps into a set of dominant modes with spatial functions and associated time-varying amplitudes, we employed EOF analyses on the 59 years’ time series (Figures 11B, C). The first EOF mode of diatom proportion explained 32.4% of the total variance. The spatial coefficients are positive over most of the inner shelf, with the most significant spatial variance in the Jiangsu coastal. There Changjiang Estuary also experienced remarkable changes. The temporal variability is shaped as an increasing trend from the end of the 1960s to the middle of the 1980s and oscillates back and forth. There has been a 10-years decline since 2000 (Figure 11C). The results of the EOF spatial coefficient and the temporal amplitude suggest that the diatom proportion change could be regional. Some areas further offshore on the outer shelf with negative spatial coefficients (Figure 11B, light blue color) suggest a diatom decreasing trend along the time series. This is because the significant changes in riverine nutrients have not affected these areas; other factors, i.e., temperature increasing could result in such variability.

FIGURE 11

Figure 11 (A) Model simulated diatom proportion climatology calculated based on the model output from 1960 to 2018. (B) The spatial coefficient of the first empirical orthogonal function (EOF) mode of diatom proportion. (C) The temporal amplitude of the first EOF mode.

Discussion

Based on a delicate grid design that can refine different spatial scales from the river channel to the estuary and the shelf, and considering various biogeochemical processes influencing the phytoplankton dynamics, our model was able to reproduce the multi-decadal changes in the nutrient distribution and the associated phytoplankton variability in the river-estuary-shelf continuum of the Changjiang Basin. In the river channel and the river mouth, the light was the main factor limiting the uptake of nutrients by phytoplankton due to a high sediment load. Thus, the nutrient concentration distribution in this area was mainly affected by riverine nutrients. Our model successfully showed that the time series of nutrients at site 1 located at the river channel were close to the riverine nutrient observations at Xuliujing station in Zhang et al. (2021). The modeled nitrogen load increased five-fold from 1970 to 2015. Regarding the variation in silica, our model captured its stable variability at a concentration between 100 and 120 μM at the river channel, which is in the same range of variance as the measurement in Zhang et al. (2021) during the same period. For the sites over the estuary, our model results suggest that nutrient dynamics are influenced by the dilution due to abundant river discharge (Wang et al., 2014) and the phytoplankton uptake. The simulated nutrient decrease from the river mouth to the inner shelf was consistent with observations conducted at a similar location (Zhang et al., 2007; Jin et al., 2010; Liu et al., 2014). Our model can also present the nutrient element variabilities in different regions. As shown in Figure 4, the phosphate in the inner shelf (site 3) did not show a prominent increasing trend compared with the other riverine sites (i.e., the river channel and the river mouth). This was different from the nitrate variability that all three sites indicated a significant increase in nitrate. The difference is that nutrient at site 3 is much more influenced by phytoplankton activity than sites 1 and 2, making the nutrient concentration processes more complicated. Moreover, the Changjiang Estuary is phosphorus-limited (Li et al., 2007; Duan et al., 2008; Ge et al., 2020a), and the increasing riverine phosphate will be consumed by phytoplankton in response, thereby slowing the rate of increase of phosphate at site 3. The nutrient change map confirmed that a significant increase in phosphate only existed in the river channel within the maximum turbidity zone; however, the nitrate increase could reach further offshore because of abundant leftover nitrate.

Based on the reproduced nutrient time series, our model simulated the community structures of the primary producers in the estuary. A significant increase in primary production during the 59 years was presented. Our model also successfully reproduced the change in the phytoplankton functional type during this period. The diatom proportion is relatively low In the river channel and river mouth due to light limitation compared with the area beyond the estuarine turbidity maximum zone. The EOF analysis suggested a significant increase in the diatom proportion along the 59-year time series. It is worth noting that there was a period after 2000 in which the diatom proportion decreased (Figures 7A, 11C). This could be associated with a decline in riverine silica because of decreases in the sediment load after the 2000s. Such asymmetric changes in nutrient stoichiometry could result in a gradual shift to silicate limitation. Our modeled phytoplankton community generally shows a similar spatial distribution as described in Liu et al. (2016), with higher diatom proportion in the estuary and inner shelf. For the temporal variability, our modeled diatom proportion decrease are consistent with the findings in Jiang et al. (2010) and Xiao et al. (2018), it reveals that there are more dinoflagellate blooms in the Changjiang Estuary, especially after 2011.

Multivariate indicators of Changjiang Estuary ecosystem

Although there are various hypotheses about the mechanisms underlying the dramatic ecosystem changes in the Changjiang Estuary, the anthropogenic impact is believed to be the main cause. Agriculture and industry developed rapidly in the Changjiang Basin since the 1980s, resulting in a 30-fold increase in nitrogen fertilizer use from 1958 to 1985, a 10-fold increase in phosphorus fertilizer from 1970 to 2010, and a drastic increase in wastewater (Duan et al., 2007; Powers et al., 2015). The amount of fertilizer use increased while the utilization efficiency continued to decrease; this resulted in nitrogen and phosphorus accumulating in the soil and rising their loads into the waters (Bouwman et al., 2013; Liu et al., 2015; Yuan et al., 2018). Due to industrial wastewater treatment, nitrogen and phosphorus from discharged industrial sewage account for only 1/3 of agricultural wastewater (Bouwman et al., 2005). Therefore, fertilizer use is the main reason for the increasing nitrate and phosphate concentrations in Changjiang and the adjacent coastal areas (Jiang et al., 2010). Before the 1970s, the concentration of nutrients in Changjiang remained low. Since the 1980s, the nitrate and phosphate levels increased rapidly (Duan et al., 2007; Powers et al., 2015). Our model could reproduce nutrient variabilities associated with these activities, e.g., the nitrate and phosphate concentrations increased by 5- and 10-fold, respectively, from 1960 to 2018, with an evident increase starting in the 1980s (Figures 3A, 4A).

Dam construction is also a cause that should be mentioned. Since the 1950s, the increased damming has led to sediment retention and reduced the suspended sediment concentration in Changjiang (Yang et al., 2011; Yang et al., 2018a; Liu et al., 2020); (Guo et al., 2020). According to relevant research, the annual sediment load at the Yichang station decreased by 97% from the 1950s to the 2010s (Guo et al., 2020), and at the DT station, it decreased by 16% from the 1960s to the 2010s (Figure 2A). Due to hydrological engineering, especially the impoundment of the TGD, the annual sediment load decreased, which led to reduced turbidity, decreased light limitation, and thus increased primary production in the Changjiang Estuary (Wang et al., 2017). In addition to affecting the physical conditions in the river and the estuary, the decline in sediment discharge also causes silicate retention (Li et al., 2007; Ding et al., 2019). With increasing nitrate and phosphate and constant silicate, the ratio of nutrients changed, which caused phosphorus and silicate limitations. The dominant phytoplankton species shifted to dinoflagellates that do not need silicon for growth (Zhou et al., 2008; Jiang et al., 2014).

A global warming-induced sea surface temperature (SST) increase is another essential factor that should be considered. Temperature changes and eutrophication are known to affect phytoplankton communities, and when simultaneous changes occur in both temperature and nutrient load, the effects become more complicated. A broadly accepted conclusion from previous studies is that the marginal seas of China experienced a robust, persistent SST increase in the last few decades until the late 1990s (Xie et al., 2002; Wang and Wang, 2013; Pei et al., 2017), followed by a significant cooling trend after, which was associated with a so-called regime shift (Kim et al., 2018). Diatoms prefer lower temperature and higher nutrient concentrations, whereas dinoflagellates are less sensitive to temperature and nutrient variations but tend to prevail under low phosphorus and high N:P ratio conditions (Letelier et al., 1993; Xiao et al., 2018). Previous studies have reported that the phytoplankton in the Changjiang Estuary have shifted to subtropical-tropical and high-salinity species to some degree, e.g., Temoraturbinata, a tropical-subtropical, high-salinity copepod species, became the dominant species in the summer of 2004 (Jiang et al., 2014). Such species adapted to a warm environment are also spreading (Jiang et al., 2017; Jiang et al., 2018b), and species that could not adapt to high temperatures migrated northward. Based on the EOF analysis of the modeled diatom proportion, we saw the most significant increasing trend was at the coast of Jiangsu, which suggests that in addition to the role of nutrients, lower temperature lower temperatures may also be responsible for the diatom increase. This variability is very close to the temporal variability of silicon (Figure 5) and SST (Kim et al., 2018). However, we couldn’t identify which was the dominant factor for diatom variability in this study. More diagnostic runs are needed in future studies.

Comparison with other estuaries

Due to continued intense anthropogenic impacts worldwide and global climate change, giant rivers and their corresponding estuaries have changed significantly since the 1950s. Eutrophication is a widespread phenomenon. The concentrations of N and P in the Mississippi River have increased drastically since the 1950s (Lohrenz et al., 2008; Smits et al., 2019). Silicate has decreased due to the burial of freshwater diatoms since 1960 (Rabalais et al., 1996). Since the 1950s, nutrients in the Chesapeake Bay have increased and resulted shifts of nutrient structure and phytoplankton species (Prasad et al., 2010; Harding et al., 2015). The Black Sea also experienced eutrophication; however, after the dam construction in the 1970s, the silicate load decreased by over 60%. The influence of Si/N ratio on phytoplankton structure was more considerable than that of eutrophication. The phytoplankton community shifted from diatoms to coccolithophores and flagellates (Humborg et al., 1997; Yunev et al., 2007; Daskalov et al., 2017). For the Nile River, the construction of the Aswan High Dam in 1964 caused a significant decline in sediment, discharge, and silicate flux. However, the Nile River did not experience eutrophication at that time because of limited nutrient inputs (Chen et al., 2021). Therefore, phytoplankton and fish populations continued to decline because of nutrient limitation. The ecosystem recovered until the early 2000s because of increasing nutrient input (Chen et al., 2021). All these studies indicate that most giant rivers and their corresponding estuaries have suffered eutrophication since the 1950s and the 1960s because of increasing N and P inputs from anthropogenic activities. The construction of dams caused the concentration of silicate to decrease, resulting in a shift in the structure of the nutrients and phytoplankton communities.

Uncertainties

We used the FVCOM-ERSEM biogeochemical model to simulate long-term nutrients and phytoplankton variability in the Changjiang Estuary. Some limitations of our study require further improvement. First, besides the simplified the biogeochemical processes described by the model, the model boundary inputs of nutrients were also simplified due to limited data sources. For example, the river discharge was based on the historical data at DT station which is located upstream from the river mouth. The riverine nutrient flux could be biased because we did not consider the river discharge and nutrient inputs along the drainage basin from DT to Shanghai; thus, the critical roles of megacities located at Changjiang deltas were not properly considered. Although our model can simulate the multi decadal variabilities, its ability to represent the interannual anomalies is limited. As we know, strong El Niño events are usually followed by massive summer monsoon over the Changjiang river basin, e.g., the crucial summer water catchments were followed by the peak of two major El Niño events in 1982/1983 and 1997/1998. These two massive summer flooding events were clearly seen from the annual water runoff time series (Figure 2) which has been used as the river boundary input in the model. However, the relationships between our simulated nutrient and phytoplankton time series and ENSO signal are less obvious. Although the El Niño signal can be seen in the physical forcing, the phytoplankton growth is regulated by both the physical forcing and nutrients in the model, thus, lack of riverine nutrient observational data in these years (e.g., 1982/1983) could lead to the bias of our simulated phytoplankton variabilities. Considering the model limitation and the main focus of this study is on multi-decadal variability, we didn’t claim the ENSO influence on the ecosystem variability.

In this study, we did not consider the long-term variability of nutrient inputs from other open boundaries, e.g., the Kuroshio, which flows on the external edge of ECS and transports a massive amount of warm, saline water and nutrients into the ECS shelf. Regarding the underlying mechanisms for long-term ecosystem change, more model sensitivity runs are needed to specify how these mechanisms act individually or in concert in this interlinked large River Delta. We also need more observational datasets to verify our simulated results and examine other important aspects of the Changjiang Estuary ecosystem with our model, such as hypoxia, acidification, and zooplankton community status.

Summary

Using a physics-biogeochemistry coupled FVCOM-ERSEM model, we simulated 59 years of ecosystem variations in the Changjiang Estuary. Our model was driven by long-term Changjiang nutrient observations and realistic atmospheric forcing and coupled with key hydrodynamics processes in the river. The model successfully simulated the nutrient concentration gradient from the river channel to the inner shelf and described how these variabilities lead to the temporal and spatial variabilities of nutrients and associated phytoplankton distributions on the adjacent shelf. Our model simulations further supported the results of nutrient ratio variabilities and associated shifts of phytoplankton species in previous studies in the Changjiang Estuary, which are based on discontinuous observations with limited temporal and spatial scales. Our coupled hydrodynamic-biogeochemical model could compensate for the deficiency of limited observations, especially in the study of long-term spatial and temporal variabilities in the ecosystem.

Data availability statement

Publicly available datasets were analyzed in this study. This data can be found here: http://xy.cjh.com.cn/.

Author contributions

SS: Methodology, Formal analysis, Investigation, Writing Original Draft. YX: Conceptualization, Formal analysis, Investigation, Writing-Review and Editing. WL: Data curation, Visualization. JG: Supervision, Project administration, Writing-Reviewing and Editing, Funding acquisition. All authors contributed to the article and approved the submitted version.

Funding

This research is supported by the National Natural Science Foundation of China (Grant No. 41606025,41776104), the National Key R&D Program of China (Grant No. 2016YFA0600903).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Bennett E. M., Carpenter S. R., Caraco N. F. (2001). Human impact on erodable phosphorus and eutrophication: A global perspective. BioScience 51 (3), 227–234. doi: 10.1641/0006-3568(2001)051[0227:Hioepa]2.0.Co;2