1 INTRODUCTION

In recent decades, especially after China's reform and opening up, the rapid rise of urban agglomerations in the Changjiang River Basin has put tremendous pressure on the environment of the Changjiang River estuary (CRE) and its adjacent waters. Large amounts of domestic sewage and industrial and agricultural pollution are discharged into the CRE and its adjacent areas, making the region one of the most eutrophic in China (Zhang, 2002). These pollutants have resulted in an increase of inorganic nitrogen and active phosphate in the area. The average annual dissolved inorganic nitrogen (DIN) flux into the sea increased from 157.0×10³ t in 1962–1977 to 900.0×10³ t in 1988– 1990, and the average annual flux of dissolved inorganic phosphorus (DIP) increased from 2.8×10³ t in 1966–1977 to 5.8×10³ t in 1978–1984 (Zhang et al., 2003). Moreover, the nutrient components (N, P, Si) and nutrient ratios in the CRE and its adjacent areas also dramatically changed due to pollutant inputs. The N꞉P ratio increased, whereas the Si꞉P ratio and Si꞉N ratio decreased significantly (Zhou et al., 2008; Chai et al., 2009). In the CRE and its adjacent areas, N is the major component for growth of phytoplankton, P is a potential limiting factor for primary production, Si is another important macronutrient, especially for the growth of diatoms, and will not become a limiting factor for phytoplankton growth in the near future although a decreasing trend of silicate might occur after the construction of Three Gorges Dam (Li et al., 2007; Duan et al., 2008; Zhou et al., 2008; Chai et al., 2009). Eutrophication and changes in nutrient components induce adverse effects, such as harmful algae blooms (HABs), depleted dissolved oxygen, and the loss of submerged aquatic vegetation and benthic fauna, thereby exacerbating ecological degradation and affecting seawater quality, ecosystem health, and human use of resources (Zhu et al., 2011; Wang et al., 2018). Therefore, determining the formation, change trajectories and ecological effects of eutrophication in the CRE and its adjacent waters is of practical significance for follow-up management and development (Hughes et al., 2005). However, most data we obtained from field investigations of the CRE and its adjacent waters starting from the mid-1980s or even later, some survey methods, research areas, and survey cycles have been inconsistent (CWRB, 1999; Wang, 2006a; Xin, 2014; Zhu et al., 2014; Wang et al., 2018). Therefore, the eutrophication history in the CRE and its adjacent area over a longer time scale remains unclear. Identifying how present-day conditions depart from the past and setting proper remediation policies to achieve an ideal ecological status in the absence of paleoecological data are difficult (Whitlock et al., 2010). Moreover, the lack of data also restricts the effective management and protection of resources in the CRE and adjacent sea areas. In this case, palaeoecological methods, which rely on proxies preserved in the sediment stratigraphy, can be a tool for providing information about the causes of change in ecosystem structure and function such as those from long-term nutrient enrichment (Cooper et al., 2004). These archives can be opened by analysing the sediment itself (its physical and chemical structures) and the remains of organisms preserved in the sediment, e. g., diatoms, dinoflagellate cysts, cladocerans, pollen, and plant macrofossils (Battarbee et al., 2001; Bennion and Battarbee, 2007; Andrén et al., 2017). By dating the sediment using, for example, radioactive isotopes (²¹⁰Pb, ¹³⁷Cs, 14C), the information contained in the sediment can be transformed into a historiography of a water body (Appleby, 2001).

Diatoms are important primary producers in the CRE and its adjacent waters. Previous studies have shown that diatoms are sensitive to variations in the environment and can be used to indicate the intensity of human activity and climate change in aquatic ecosystems (Battarbee et al., 2001; Andrén et al., 2017). Moreover, diatom cell walls deposited in sediments can be preserved well enough to use for species identification and numeration (Reavie et al., 1995; Rosén et al., 2000). Therefore, palaeoecological techniques based on diatoms can be used to quantify and infer trends in the trophic status in the CRE and its adjacent areas over long periods of time. Although palaeoecological techniques and transfer functions were originally developed for and applied to lake (Reavie et al., 1995; Rosén et al., 2000; Zou et al., 2021) rather than estuarine environments, such methods have been proven to be applicable to coastal environments, including estuaries (Andrén, 1999; Clarke et al., 2003; Weckström et al., 2004; Clarke et al., 2006; Ellegaard et al., 2006; Weckström, 2006; Logan et al., 2011; Andrén et al., 2017). A study establishing diatom-based transfer functions for retrospective assessment of water eutrophication conditions in the CRE and its adjacent areas has not yet been reported.

Due to the complex hydrodynamic conditions in the CRE and its adjacent areas, both physical and chemical elements change significantly with seasons and regions, resulting in differences in diatom community structures (Song et al., 2008; Liu et al., 2015; Wang et al., 2016). Contemporary relationships between the distribution of diatom assemblages and environmental variables in the CRE and adjacent areas suggest that DIN is one of two important factors (Fan et al., 2019). Variations in diatom assemblages stored in stratigraphic records may reflect historical changes in DIN in the CRE and its adjacent areas; this finding could be the basis for developing an index to obtain an independent record of past changes in DIN in this region.

In this study, we developed a diatom-based DIN record for the past 50 years from a sediment core (core DH8-2) taken from the Fujian-Zhejiang coastal mud area south of the CRE. This study had two main purposes: ⅰ) to quantitatively reconstruct the DIN changes in the southern CRE and ⅱ) to explain the variation in nutrients and the potential influence of human activities in the CRE and its adjacent waters by combining modelled inferences of DIN and biogenic elements.

2 MATERIAL AND METHOD 2.1 Study area

The CRE, one of the largest estuaries in the world, is influenced by complex hydrological conditions, including the Changjiang Diluted Water (CDW), the Taiwan Warm Current (TWWC), continental river runoff, and the Kuroshio branch current in the East China Sea (Fig. 1) (Ning et al., 2011; Yang et al., 2012). Therefore, the distribution and variation in temperature and salinity in the CRE and its adjacent areas exhibit seasonal and regional differences. In summer, the characteristics of CDW are prominent. The low-salinity surface water can reach 123°E and has a trend of east expansion off the CRE, which results in high temperature and low salinity in the estuary area, and low temperature and high salinity in the area outside of the estuary. In winter, low temperature and low salinity are present in the estuary, and high temperature and high salinity exist in the area outside the estuary. Compared with the distributions of salinity in summer, the salinity isolines move inside, and the low salinity water tongue is obviously shrinking (Tang and Wang, 2004; Chen, 2009). The circulation system also provides power for the transportation of foreign materials in marine sediments and establishes a suitable living environment and nutrient supply for organisms. The sediment types in the CRE and its adjacent areas are mainly silty clay, and the content of sand particles is generally low within 10% of the coastal waters and does not change significantly (Xu et al., 2016). Two major mud areas, the near-shore mud area (including the CRE mud area and the Fujian-Zhejiang coastal mud area) and the south-western mud area of distal Jeju Island, are formed by the deposition of suspended matter carried by the Changjiang River (Guo et al., 2000). In the Fujian-Zhejiang coastal mud area, the sources of sediment are mainly the Changjiang River, Qiantang River, and Oujiang River, etc. (Yang et al., 2016). In the past 100 years, the sedimentation rate in the Fujian-Zhejiang coastal mud area of the East China Sea ranges from 0.79–3.34 cm/a, on average of 1.97 cm/a, making it a sub-high sedimentation intensity area (Shi et al., 2010). Due to its stable sedimentary environment and high deposition rate, sediments could truly and effectively reflect the changes in the paleoenvironment and the impacts of human activities in this area (Guo et al., 2000; Xiang et al., 2006).

2.2 Core sampling and pretreatment

A 50-cm long sediment core named DH8-2 was sampled from the Fujian-Zhejiang coastal mud area (27.2°N, 121.5°E; Fig. 1a) at water depth of 58 m during the 2012 spring open cruise of the Institute of Oceanology, Chinese Academy of Sciences (IOCAS). The core was generally greyish brown, and based on the observation of lithology, the lithostratigraphy in core DH8-2 was homogeneous, and mainly comprised silt and clay. A small number of shell fragments were found at 20–25 cm and 30 cm. After cutting intervals every 1 cm, subsamples were sealed in polyvinyl chloride (PVC) pouches, stored at -20 ℃, and then brought back to the laboratory for further analysis. To determine the total organic carbon (TOC), total nitrogen (TN), and biogenic silica (BSi), two parallel samples were taken every 1 cm, for a total of 21 samples. Samples for diatom and dating determination were collected every 4 cm, for a total of 11 samples, as shown in Table 1.

2.3 Diatom

The sediment samples were processed for diatom analysis as per Battarbee et al. (2001). Approximately 10 g of freeze-dried sediment samples were treated with 10% HCl, 30% H₂O₂, and ammonium hydroxide to completely remove carbonates, metal salts, metal oxides, and organic matter, and to reduce the impact of clay. Then, siliceous microfossils were suspended by adding zinc bromide (specific gravity: 2.4) with centrifugation at 2 000 r/min for 10 min. The suspended valves were diluted to certain volumes, and aliquots were transferred to cover slides. After air-drying of the aliquots, permanent slides were prepared with Naphrax^TM (1–2 drops). The diatom valves were counted in three slides prepared from each sample. Slides were examined under an Olympus Ⅸ71, and at least 300 valves were counted in each sample and identified according to Hasle and Syvertsen (1996) and Guo and Qian (2003). Some diatom taxa were identified only to the genus level instead of the species level due to slide and identification issues. One silicoflagellate species (Dictyocha fibula) was identified, but due to its low relative abundances, it did not appear in the subsequent analyses.

The diatom absolute abundance (D_abs) was expressed as the number of valves per gram of dry sediment (valves/g DW) and was calculated as follows:

where N is the total number of diatoms in three slides; V is the total volume of diluted suspended diatoms; V₁ is the volume of diluted suspended diatoms applied on three slides and W is the weight of the dry sediment sample.

The species diversity and richness indexes were calculated according to the Shannon-Weaver index (Hʹ) (Shannon and Weaver, 1949; Margalef, 1968).

where Hʹ, Shannon-Weaver diversity index; D, species richness index; S, number of species; N, total individuals of all species; P_i, ratio of the individuals of species i to the total individuals of all species.

2.4 Geochemistry and dating

The determination of TOC and TN was performed on a Vario MACRO cube element analyzer according to Sun et al. (2014); the precision of the instrument was < 0.5%, and the recovery efficiency was > 99.5%. An improved wet alkaline method was employed to extract BSi using Na₂CO₃ and a sequential extraction technique reported by Wang et al. (2014), and then the dissolved silica content in the extractions was measured using the molybdate blue spectrophotometric method (Mortlock and Froelich, 1989; Conley, 1998).

The sediment samples were dated using ²¹⁰Pb and ¹³⁷Cs dating techniques, which are widely used for marine and lake sedimentary environments on a hundred-year scale (Deng et al., 2006). Measurement of ²¹⁰Pb and ¹³⁷Cs followed DeMaster (1981); a 2–55-g freeze-dried sample was ground (< 63 μm), wax-sealed in plastic vials, and equilibrated for 15 days. Next, the samples were analyzed using high-purity germanium gamma spectrometry (ORTEC). The sedimentation rate of the present study is based on a constant rate supply model of ²¹⁰Pb (DeMaster, 1981; Zaborska et al., 2007) and ¹³⁷Cs activity.

2.5 Numerical method

The calibration data set (including diatom data from 34 surface sediments and 8 water environmental variables) used for the development of the transfer function was collected in 2015 (Fan et al., 2019). The surface sediment samples were sampled using a gravity corer or grab, and the top 3 cm represented recent sediment deposition (Fig. 1a). Water samples were collected from the surface, intermediate, and bottom layers before sediment sampling at each sampling site seasonally throughout 2015, and then the mean annual values of each water layer were used for the training data set. Stations were selected to cover a relatively long nutrient gradient, with salinity ranging between 24.96 and 34.15, DIN ranging between 8.75 and 34.63 μmol/L, DIP ranging between 0.28 and 0.89 μmol/L, and dissolved inorganic silicate (DSi) ranging between 7.77 and 35.50 μmol/L. A total of 114 diatom species were identified in the modern assemblages, and as diatoms are rare within the samples and lack adequate data to characterize any environmental associations, 37 taxa accounting for > 2% of at least two samples were used in the development of the transfer function. Prior to analysis, the diatom taxa percentage data were square root transformed to decrease the weighting of dominant species (Bigler et al., 2006). Eight water variables were lg(x+1) transformed to eliminate the impact of their units. A thorough description of the methodology can be found in Fan et al. (2019).

Partial redundancy correspondence analysis (RDA) was performed in the R software environment, and revealed that DIN was independently associated with the diatom composition and distribution; therefore, taking DIN as the only environmental variable, the results of detrended canonical correspondence analysis (DCCA) showed that the gradient length of DIN was 2.73 standard deviation (SD) units, which is greater than 2.4 SD units, indicating that a unimodal response model was more appropriate for the reconstruction (Hassan et al., 2009). In addition, the ratio of the first constrained eigenvalue (λ₁) to the second unconstrained eigenvalue (λ₂) in a constrained RDA showed that DIN was the dominant gradient with a λ₁/λ₂ of 0.83, indicating the robustness of the training set for reconstructing DIN. Therefore, the transfer function was developed to reconstruct DIN concentrations using a weighted average (WA) model and weighted average partial least squares (WA-PLS) model. Weighted averaging regression is the most widely used numerical method for developing diatom-based transfer functions. The method was first used in paleolimnology by ter Braak and van Dam (1989) and developed by Birks et al. (1990) during the Surface Water Acidification Project (SWAP). It has subsequently been used for the reconstruction of a wide range of variables including salinity (Hassan et al., 2009; Li et al., 2012; Volik et al., 2017), TN (Clarke et al., 2003; Weckström, 2006; Andrén et al., 2017), TP (Enache and Prairie, 2002; Logan et al., 2011; Murnaghan et al., 2015; Witak et al., 2017; Sochuliaková et al., 2018; Sui et al., 2020), pH (Enache and Prairie, 2002; Sienkiewicz and Gąsiorowski, 2017; Narancic et al., 2021; Sienkiewicz et al., 2021), and sea level (Horton et al., 2006; van Soelen et al., 2010; Gomes et al., 2017). The performances of the models were assessed by leave-one-out cross-validation (999 permutation cycles) (Birks, 1995). The model with strong performance was identified as having a combination of a high coefficient of determination (R²) between observed and predicted values, low mean and maximum biases, and a low root mean squared error of prediction (RMSEP). Transfer functions were determined using the computer program C2 (Juggins, 2007).

A cluster analysis was performed on diatom percentages to define diatom zones using the computer program SPSS 21.

3 RESULT 3.1 Diatom-DIN transfer function and evaluation of its reliability

The WA and WA-PLS models were used to develop a diatom-DIN transfer function based on revised modern calibration data from the CRE and its adjacent areas (Fan et al., 2019). According to the performances of several transfer functions (Table 2; Fig. 2), the five-component WA-PLS model, which had the highest R_jack-knifed² value (0.822) and the lowest RMSEP_jack-knifed (0.142 lg DIN μmol/L), appears to perform better than the WA model based on the training data set with 34 samples. Therefore, WA-PLS with five components was used for reconstructing DIN concentrations.

To further evaluate the reliability and applicability of the WA-PLS diatom-DIN transfer function, we selected 21 additional sets of surface sediment diatom data and corresponding upper water DIN data (collected in July 2012, unpublished data) to verify the model. Among these data, 18 samples were collected from the CRE and its adjacent areas (Fig. 1b, rhombuses, C1-C18), and 3 samples were collected from the Yellow Sea (Fig. 1b, triangles, Y1–Y3).

The verification results show that among the 18 samples in the CRE and its adjacent areas, the agreement between the diatom-inferred DIN (DI-DIN) data and the actual measured value is close, with R² values up to 0.85 (Fig. 3), which is consistent with the results obtained by the jack-knife technique in the previous model. However, the agreements between the DI-DIN data and the actual measured value of the Yellow Sea are poor (Table 3), especially for the northernmost Y1 sample, where the inferred DIN is nearly 13 times higher than the observed value. The above results prove that the diatom-DIN quantitative transfer model based on the calibration data from the CRE and its adjacent areas is very likely unsuitable for other research areas.

3.2 Age model

Eleven ²¹⁰Pb_ex (excess ²¹⁰Pb) and ¹³⁷Cs dates were obtained from core DH8-2 (Table 1). The ²¹⁰Pb_ex activity profile in the sediment core showed an overall pattern of exponential decline (Fig. 4a), and a constant initial concentration (CIC) model (Appleby and Oldfield, 1983) was used to calculate the deposition rate of the core sediment. According to the logarithmic linear fitting of ²¹⁰Pb_ex and depth, the average deposition rate was 0.996 cm/a. Therefore, the core sediments can be speculated to contain a total of 50 years of deposition records from 1962 to 2012. The profiles of ¹³⁷Cs in the core sediment revealed two major peaks at core depths of 25–26 cm and 35–36 cm (Fig. 4b), which probably correspond to 1986 (the Chernobyl nuclear accident) and 1976 (China's nuclear testing), respectively. The ²¹⁰Pb chronology for this core was in general agreement with the ¹³⁷Cs chronology.

According to the CIC model and ¹³⁷Cs chronology, the sedimentation rate was relatively stable and fluctuated little, remaining at approximately 10.0 mm/a. This may be because core DH8-2 was collected in the central and southern parts of the Fujian-Zhejiang coastal mud area, south of the CRE. The hydrodynamic conditions and terrestrial matter input of the Changjiang River, the Minjiang River and other coastal rivers are relatively stable.

3.3 Geochemistry

The vertical distributions of TOC, TN, and BSi were roughly the same, showing a general upwards trend from bottom to top (Fig. 5). The values of the three indexes were stable and relatively low from the bottom layer to 35 cm, above which the variation in the three indexes showed fluctuations and relatively high values (Fig. 5). From the bottom to 35 cm, the average contents of TOC, TN, and BSi were 0.230%, 0.023%, and 0.143%, respectively, which were the lowest in each segment. From 35 cm to the top of the core, the TOC, TN, and BSi contents significantly increased, with high average contents of 0.408%, 0.038%, and 0.414%, respectively. The correlation analysis showed a significant correlation between the BSi and TOC contents (P < 0.01); the BSi/TOC ratios varied from 0.368 to 1.481 with an average of 0.909 (Fig. 5), much higher than the Redfield ratio (BSi/TOC=0.151), indicating that diatoms contributed significantly to total primary productivity in this region (Tréguer et al., 1995; Lü et al., 2010). The range of TOC/TN ratios in the core profile was 7.37–20.67, with an average of 11.17. This result suggests that the sources of core DH8-2 were mainly terrestrial organic matter and that marine materials played a relatively minor role in this area. Core DH8-2 was collected in the central and southern parts of the Fujian-Zhejiang coastal mud area, south of the CRE, and the sediment in this area was mainly (exceeding 80%) transported from the Changjiang River to the sea by the Fujian-Zhejiang coastal current (Xiao et al., 2005). The Oujiang River and Minjiang River also contributed sediment to this area to a certain extent.

3.4 Diatom

Diatom assemblages were enumerated from samples taken from 11 depths and were abundant and well preserved. A total of 70 taxa, consisting of 38 centric diatoms and 32 pennate diatoms, were identified. The diatom assemblages were mainly marine taxa and a few brackish-freshwater and

freshwater taxa. Diatoms for planktonic and benthic life forms can be seen throughout the core. These diatom taxa belonged to 17 families and 30 genera, of which Coscinodiscus (16 species), Nitzschia (6 species), and Diploneis (4 species) were the main genera with most diatom taxa. Diatoms with an abundance of more than 2% in at least one stratigraphic level are shown in Fig. 6. Some clear trends in diatom stratigraphy are identifiable as shifts in life forms (planktonic and benthic) and nutrient requirements. The core was divided into 3 diatom assemblage zones (Zones Ⅰ, Ⅱ, and Ⅲ) based on the cluster analysis and described as follows:

Zone Ⅰ (50–40 cm, ~1962–1972) represents the bottom part of the core and was characterised by the dominance of benthic diatoms, e.g., Podosira stelliger, which accounted for more than 50% of the diatom abundance, and the frequencies of tychoplanktonic taxa, e.g., Coscinodiscus argus and Cyclotella striata were also high, with abundances of approximately 10% and 20%, respectively (Fig. 6). This zone had the lowest diatom absolute abundance (Dabs), Shannon-Weaver index (H′), and richness index (D), with averages of 1.42×10³ valves/g DW, 2.12, and 1.73, respectively (Fig. 7), indicating low diatom community diversity in this period.

Zone Ⅱ(40–10 cm, ~1972–2004) is characterized by a distinct increase in the frequency of planktonic taxa and dramatic decreases in the abundances of P. stelliger and C. striata, although both species still represented a numerically important part of the diatom community, with average frequencies of 14.24% and 5.55%, respectively. C. argus also became notably less abundant. Paralia sulcata and a markedly high abundance of planktonic taxa, such as Cyclotella stylorum (0–26%) and Actinocyclus ehrenbergii (0–12%), also occurred in this zone (Fig. 6). The species diversity in this zone increased significantly compared with that in zone Ⅰ, which was the highest among the three zones, and the average Shannon-Weaver index (H′) and richness index (D) reached 4.29 and 3.72, respectively (Fig. 7).

Zone Ⅲ (10–0 cm, ~2004–2012) is characterized by markedly elevated abundances of P. sulcata, Actinoptychus senarius, and Thalassiosira spp. and distinct decreases in and even the disappearance of P. stelliger, C. striata, and C. stylorum. The frequency of planktonic taxa further increased in this zone, up to 65%. The species diversity in this zone was lower than that in zone Ⅱ, but it was still high, as confirmed by the high Shannon-Weaver index (H′) (3.09–4.20) and richness index (D) (2.64–3.37) (Fig. 7).

3.5 DI-DIN reconstruction

DI-DIN concentrations were calculated based on the DIN transfer function developed by the five-component WA-PLS model. The DI-DIN (Fig. 6) indicated that the DIN concentration showed an apparent upwards trend from the bottom to the top of the core sediment. The DIN fluctuated between 5.74 µmol/L and 6.20 µmol/L in the lowest part of the core and then started to increase significantly between 1972 and 1988 and reached between 9.22 µmol/L and 10.41 µmol/L between 1998 and 2012. The maximum value (10.41 µmol/L) was reached in 2008. As shown in Fig. 6, the change trend of DI-DIN was similar to that of the abundance of planktonic/benthic species in the diatom community, and the correlation analysis showed a significant positive correlation between the DI-DIN values and the ratio of planktonic/benthic species (R²=0.92, P < 0.01); that is, with the increase in DIN in water, the number of planktonic species in the diatom community increased, and the abundance of benthic species decreased.

4 DISCUSSION 4.1 Application of the diatom-DIN transfer function in the CRE and its adjacent areas

In this study, the WA-PLS model performed better than the WA model in developing the DI-DIN transfer function. The WA model ignored correlations in the diatom species data that remained after fitting the DIN variable; the diatom species response in the WA model was also influenced by secondary or other environmental variables, such as salinity and temperature, that were not taken into account, while in the WA-PLS model, successive components were extracted from the training data set to increase the predictive power of the regression model (ter Braak and Juggins, 1993; Birks, 1995). For the DI-DIN transfer function based on the 5-component WA-PLS, the R_jack-knifed² value was 0.822, and the RMSEP_jack-knifed value was 0.142 lg DIN μmol/L, indicating the powerful predictive ability of this function (Table 2). Two possible explanations exist for this finding: 1) The 34 sites contained in the training data set spanned a wide range of DIN concentrations, between 8.75 and 34.63 μmol/L (Fan et al., 2019), and the diatom composition had a significant fit to a symmetrical unimodal distribution along the DIN gradient, therefore yielding a high correlation between the observed and inferred nutrient values (Potapova et al., 2004). 2) A low frequency of benthic diatoms was included in the DI-DIN transfer function. Jones and Juggins (1995) showed that many benthic diatom species had no direct response to the nutrient concentration in the upper water body, which was mainly affected by other environmental variables, such as water depth, lake transparency, and macrophyte cover/biomass; therefore, the low frequency of benthic diatoms may have improved the effectiveness of the DI-DIN transfer function to a certain extent (Sayer, 2001).

Despite this good performance, the model should be interpreted with caution because the accuracy of DIN reconstructions and their paleoecological interpretation may be affected by other factors. First, the data set was relatively small (Fan et al., 2019), and other environmental variables not included in the training data set may be important in the distribution and composition of the diatom community. Therefore, the agreement between fossil and modern data sets could be improved by incorporating more environments into the model. Second, the distributions of diatom assemblages in estuaries and coastal areas are influenced by taphonomic alterations related to physical mixing and currents (Hassan et al., 2009). However, Hassan et al. (2008) showed that diatom assemblages reflect major estuarine and coastal environmental gradients, although mixing of allochthonous and autochthonous diatoms was shown to occur. Many studies have shown that diatom transfer functions appear to overestimate or underestimate the variable being modelled in reconstructions at measured values that are too high or too low (Schönfelder et al., 2002; Weckström et al., 2004, 2007; Dong et al., 2006; Weckström, 2006; Andrén et al., 2017). Our model performed well when the DIN values were below 20.27 μmol/L (1.33 lg DIN μmol/L), and no consistent underestimation or overestimation was found (Fig. 2). This was because most of the sites in the data set had relatively low DIN concentrations, and the corresponding diatom species taxa were abundant. Comparing the DI-DIN values with the measured DIN, underestimation clearly occurred at the high DIN values (Fig. 2); this was due to the lack of sites with high DIN concentrations and indicated that the modelled results may have relatively low reliability (Bradshaw and Anderson, 2001). A diatom transfer function study from the Gulf of Finland similarly underestimated nitrogen concentrations during the most eutrophic phase due to a lack of modern analogues in the training data set (Weckström et al., 2004; Weckström, 2006). A diatom-based TN reconstruction applied in coastal waters of the Bothnian Sea also underestimated concentrations at high TN values for the same reason (Andrén et al., 2017). In general, an insufficiently high DIN concentration causes a decrease in analogue quality. Birks et al. (1990) showed that when environmental gradients become shorter, the responses of diatom species to the environment will be truncated, leading to an increase in model derivation error. Therefore, this DI-DIN transfer function must be applied in a reasonable DIN concentration range. In this study, the disadvantage of the transfer function when reconstructing high DIN concentrations did not affect our aim to quantify the low background nitrogen conditions in the CRE and its adjacent areas in recent years.

In addition, we found that the DI-DIN transfer function did not perform well in the 3 samples from the Yellow Sea area (Table 3). The possible reasons for this finding are as follows: 1) three samples could not reflect the actual verification effect of this DI-DIN transfer model in the Yellow Sea area; 2) due to differences in environmental factors, the compositions and structures of diatoms preserved in the surface sediments of the Yellow Sea were significantly different from those of the CRE and its adjacent areas (i. e., the data we used to develop the DI-DIN transfer model), as reflected in a marked decrease in the quantity and variety of Coscinodiscus sp. (Wang et al., 2016; Fan et al., 2019).

4.2 Implications of the reconstructed DIN and diatom changes

The DI-DIN profile showed a distinct increase (by more than 80%) in the DIN concentration in the Fujian-Zhejiang coastal mud area south of the CRE from 1962 to 2012 (Fig. 6), suggesting that eutrophication was intensified in this area over this time period. The terrestrial inputs in this area were mainly from the Changjiang River, Oujiang River, and Minjiang River; among them, water from the Changjiang River had a major contribution to the DIN and diatom composition, while water from the other rivers had less impact (Xiao et al., 2005). Wang (2006b) found that the DIN of water from the Changjiang River increased exponentially and by a factor of five from the 1960s (20 μmol/L) to the end of the 1990s (120 μmol/L), and a 2-fold increase in DIN concentration was observed in the surface water of the Changjiang plume. The concentrations of TOC, TN, and BSi in soils and sediments can also serve as important paleoenvironmental proxies of the primary production of phytoplankton and the eutrophication level in upper water layers (Meyers, 2003; Tallberg et al., 2015). In this study, the contents of TOC, TN, and BSi showed an overall upwards trend, from 1962–2012 (Fig. 5), consistent with the intensified eutrophication during this period. According to the DIN reconstruction (Fig. 6), changes in the DIN concentration can be divided into three stages, namely, low-level (1962–1972), rapid (1972–2004), and continuous (2004–2012), which are closely related to China's development pattern and similar to the changes in riverine concentration and flux of DIN from the Changjiang River (discussed below) (Wang, 2006a; Xin, 2014; Wang et al., 2018). Therefore, the DIN reconstruction in this study was reliable to a certain extent and can reflect the historic eutrophication changes in the CRE and its adjacent areas. Anderson and Vos (1992) and Andrén (1999) noted that increases in diatom biomass were accompanied by increases in the levels and changes in the components of nutrients, while the percentage of plankton increased and the percentage of benthos decreased. Comparable results were obtained in this study and a significant positive relationship was found between DIN and the planktonic/benthic ratio.

In the stage with low levels of change in the DIN concentration (1962–1972), the DIN was stable at approximately 6 µmol/L, suggesting a relatively low nutrient level. China's industry and agriculture were in their infancy in this period, and the Changjiang River and coastal waters were less affected by anthropogenic activities. Although the riverine concentration and flux of DIN from the Changjiang River increased beginning in the 1960s, the increasing trend was not significant in the following 10 years (Wang, 2006a). The contents of TOC, TN, and BSi were smaller before 1976 (Fig. 5), indicating relatively low primary production in water bodies with less human influence (Meyers, 2003). The low-nutrient condition was also confirmed by the low and stable numerical values of absolute abundance and diversity index of diatom fossils (Fig. 7).

In the rapid increase (1972–2004) stage, the average DIN increased nearly 40% from the last stage, ranging from 7 µmol/L to 9 µmol/L. The inferred DIN values were close to the observed values. Wang (2006b) showed that the DIN of this area was between 8 µmol/L and 9 µmol/L during 1997–1999, which further supports the reliability of the DI-DIN transfer function applied in this research. The rapidly increasing DIN suggested a relatively high nutrient concentration and a rapid worsening of eutrophication; the concentrations of three geochemical parameters, especially TN, increased significantly during the period from 1976 to the 1990s (Fig. 5), which also supports these conditions. This stage corresponded to the most rapid development period of China. In the 1970s, China started to enact policies of reform and open up to the outside world; since then, China has witnessed great progress in economic and social development, but severe environmental problems have followed. Various anthropogenic activities (e.g., industrial and agricultural production, urban living, and aquaculture) directly or indirectly result in nutrient inputs to coastal waters (Qu and Kroeze, 2010; Wang et al., 2018). The use of chemical fertilizer in China increased nearly 12-fold from 1970 (3.51×10⁶ t) to 2000 (41.4×10⁶ t), resulting in a sharp increase in the export of dissolved and particulate nitrogen from rivers to coastal waters (Qu and Kroeze, 2010; Luan et al., 2013). Sewage discharge in the Changjiang River Basin increased by 63% from the 1980s to 1999, the riverine concentration and flux of DIN from the Changjiang River entered a stage of rapid increase from the 1980s to the 2000s, and the DIN concentration in the CRE tripled from the 1980s to the 2000s (CWRB, 1999; Wang, 2006a; Xin, 2014; Wang et al., 2018). Moreover, the Si/N ratio decreased from approximately 4 to 1 over the same period (Wang et al., 2018). Increases in the nutrient levels and changes to more bioavailable nutrient types have increased eutrophication and phytoplankton biomass, which in turn have triggered frequent harmful algal blooms (Yu et al., 2017). In this stage, the absolute abundance and diversity index of subfossil diatoms preserved in sediments showed a marked increase, and the species composition was indicative of rich nutrients. For example, the abundances of A. ehrenbergii, P. sulcata, and Coscinodiscus radiatus gradually increased (Fig. 6), indicating a significant increase in nutrients during this period (McQuoid and Nordberg, 2003). The change in life forms can be attributed to changes in habitat availability. An increase in the abundance of planktonic taxa has frequently been reported as a response to eutrophication in freshwater and marine environments (Bennion et al., 2004; Cooper et al., 2004; Wang, 2006a; Xin, 2014). This trend is supported by diatom fossils in our study. The planktonic/benthic ratio was significantly higher in this stage than in the last stage (Fig. 6). Wang (2003, 2006b) reported that the Si/N ratio decreased from approximately 4 to 1 in the CRE and its adjacent areas over the same period, and Si limitation changed the species composition of diatoms. Under such conditions small-bodied taxa with low silica requirements will become more competitive. The taxa with large individuals and more heavily silicified, such as P. stelliger (cell diameter: 30– 85 μm), C. striata (cell diameter: 10–50 μm), and C. argus (cell diameter: 31–110 μm), decreased in percentage throughout the stage. The average percentages of the above three taxa were 14.24%, 5.55%, and 2.02% respectively, which was significantly lower than the values of the previous stage (59.42%, 15.22%, and 7.03%, respectively). The frequency of smaller individuals, such as P. sulcata (cell diameter: 16–36 μm), Thalassiosira oestrupii (cell diameter: 12–43 μm), and Thalassiosira simonsenii (cell diameter: 3–18 μm), gradually increased (Fig. 6). Similar phenomena can be seen in the Gulf of Mexico (Dortch et al., 2001), the Bohai Sea (Wang and Li, 2006), and Jiaozhou Bay (Liu et al., 2008).

In the 2004–2012 stage, the DIN was maintained at a relatively high level of approximately 10 µmol/L. In this period, the DIN flux in the East China Sea increased significantly, resulting in a continued rise in nutrients in the upper water body. Large fluctuations in the concentrations of TOC, TN, and BSi have been maintained since the 1990s, indicating that this area was continuously affected by human activities and that the level of water productivity was high. The structure of the sedimentary diatom community changed, suggesting more severe eutrophication in this period. For example, the frequency of planktonic taxa continued to increase, reaching a maximum of 65.48%, while the frequency of large benthic taxa, e.g., P. stelliger and C. striata, further decreased, with some species disappearing (Fig. 6). Small-bodied taxa, such as Thalassiosira spp., exhibited their highest abundances in this stage (Fig. 6). These fragile taxa are often not well preserved in sediments (Weckström, 2006). Their increased percentage may be a real response to higher nutrient concentrations. Moreover, dinoflagellates have gradually replaced diatoms to become the main harmful algal bloom (HAB) species in the East China Sea since 1998, and HABs caused by Prorocentrum donghaiense, Alexandrium catenella, and Karenia mikimotoi occur frequently (Zhou and Zhu, 2006; Zhang, 2011), further supporting the increase in nutrients and changes in the types of nutrients.

5 CONCLUSION

In this study, the sedimentary diatom record and geochemical data were used to reconstruct historical coastal DIN concentrations in the CRE and its adjacent areas. A training data set from 34 sites in the CRE and its adjacent areas was developed and used to create a diatom-based transfer function for DIN concentrations. The model was applied downcore (1962–2012) in the Fujian-Zhejiang coastal mud area south of the CRE. Core sediments allowed for the quantitative reconstruction of trophic changes inferred from sediments in the CRE and its adjacent areas.

The reconstructed DIN values of 5.74–0.41 µmol/L during the 1962–2012 period indicate that nutrient levels in the CRE and its adjacent areas were closely linked to China's development pattern and the environmental protection policies and measures in place. Distinct changes in diatom assemblages, DI-DIN concentrations and geochemistry indexes indicated substantial variability in the nutrient status of the CRE and its adjacent areas in the last 50 years. Before approximately 1970, the area had a stable and relatively low trophic level, with low TOC, TN, and BSi concentrations and low diatom species richness and was dominated by benthic taxa. Since the reform and opening up in the 1970s, nutrient status has increased significantly, accompanied by rapid economic and social development in China. The aggravation of eutrophication in the CRE and its adjacent areas was confirmed by the increase in BSi contents in sediments and changes in sedimentary diatom structures, with a marked increase in planktonic taxa and species favoring nutrient-rich environments. The results of this research show that nutrients have changed in the CRE, and the levels and types of nutrients over the time period associated with the sediment core could provide reference levels for subsequent resource management in this region.

6 DATA AVAILABILITY STATEMENT

All data generated and/or analyzed during this study are available upon request by contact with the corresponding author.

7 ACKNOWLEDGMENT

We thank the captain, crew, and scientific party of the R/V Kexue No. 3 Cruise for their assistance with the data acquisition.