Nutrient Dynamics and Ecosystem Metabolism of Megacity Rivers: Influence of Elevated Nutrient Concentrations in Beijing’s Waterways

Abstract

This study evaluated nutrient flux (nitrate (NO₃⁻), ammonium (NH₄⁺), phosphate (PO₄³⁻), and dissolved organic carbon (DOC) at the sediment-water interface and river ecosystem metabolism (REM) to investigate how these ecological functions vary in Beijing’s urban waterways. Three tributaries of the River Beiyun were selected. Water quality varied across the study sites as each receives a mixture of wastewater treatment plant (WWTP) effluents and tributary inflows. A chamber technique was applied where water-specific nutrient concentrations were measured at two exposure times (3 and 10 min). Under the actions of physical and biological processes, NO₃⁻ and NH₄⁺ flux was primarily controlled by equilibrium concentration and the N-cycle. However, bioabsorption appeared to regulate DOC flux. Specifically, NO₃⁻ flux ranged from −0.31 to +0.30 mg/(m²·s), NH₄⁺ was −0.01 to +0.05 mg/(m²·s), PO₄³⁻ was −0.01 to +0.01 mg/(m²·s), DOC was −0.04 to +0.13 mg/(m²·s). We applied the nighttime slope regression to estimate gross primary production (GPP) and ecosystem respiration (ER). Except in summer, net ecosystem production (GPP+ER) less than 0 indicated heterotrophic study reaches. Structural equation modelling revealed that nutrient dynamics and water temperature were the primary factors driving REM. Our study provides the needed systems-based understanding of vital ecological processes to improve in-stream management.

1. Introduction

Available nitrogen (N) and phosphorus (P) are important nutrients supporting aquatic life activities [1,2]. It is well known that excessive N and P cause eutrophication and subsequent water quality problems, such as hypoxia, which is a threat faced by urban rivers in most parts of the world [3,4,5]. Therefore, understanding the dynamics of N and P is an important component of aquatic conservation and ecology research [6,7,8]. In urban rivers, N and P are regarded as major pollution stressors and reach urban waterways via multiple sources [9,10,11]. For example, vehicle exhaust and atmospheric deposition [12,13], urban impervious surface [14,15], and wastewater treatment plant (WWTP) effluents [16,17] all can contribute to elevated N or P in urban rivers. The impairment of water quality and ecological function caused by excessive nutrients in urban rivers is the core concept of the Urban Stream Syndrome (USS) which provides a framework for evaluating responses to urbanization [4,18,19,20,21].

Dissolved organic carbon (DOC) is also an important indicator of river water quality [22]. Excessive DOC is a major source of organic pollution [23], while also being involved in important ecological processes, for example, it provides electrons for biogeochemical processes (e.g., denitrification) and serves as an energy source and respiratory substrate for many aquatic organisms [24,25]. The sources of DOC in urban rivers include both natural and anthropogenic processes. These can include riparian soil erosion, organic matter discharged by human activities, WWTP effluents, and biological residues produced as a result of photosynthetic processes of aquatic vegetation [26,27,28].

As an important part of river ecosystems, river sediments play a key role in nutrient dynamics and water quality and can serve as both a source and sink of nutrients [16,29,30]. The uptake and/or regeneration of nutrients by sediments depend on a variety of influencing factors, making river ecological protection more challenging. Nitrate (NO₃⁻) and ammonium (NH₄⁺) are the main forms of nitrogen in sediments [31,32], and their fluxes are primarily controlled by equilibrium concentration at the sediment-water interface and nitrogen form transformation within the N-cycle [29,33,34], including nitrification, denitrification, dissimilatory nitrate reduction to ammonium (DNRA), mineralization and anammox. These processes are mainly performed by microorganisms such as nitrifying bacteria and denitrifying bacteria, which are affected by a range of factors including water temperature, dissolved oxygen (DO) and associated microbial activity [35,36]. Phosphate (PO₄³⁻) is another macronutrient required by aquatic organisms [37]. PO₄³⁻ flux is related to the adsorption of sediment minerals, for example, iron insoluble compounds can adsorb a large amount of PO₄³⁻ under oxidation conditions [38,39,40,41], resulting in the uptake flux of PO₄³⁻. In the river restoration practice, alum is often added to water to promote PO₄³⁻ adsorption [42,43,44,45]. Additionally, it has been demonstrated that sediment particle size plays an important role in PO₄³⁻ flux, such as PO₄³⁻ release caused by resuspension of fine sediment particles [46,47,48,49,50]. Dissolved organic carbon (DOC) can be an important carbon source in rivers to support biological activities which can regulate fluxes. Furthermore, some studies have found that anions in DOC and PO₄³⁻ compete for adsorption sites, which can affect the rates of both fluxes [51,52].

Nutrient supply affects aquatic bioactivity, which in turn affects gross primary production (GPP) and ecosystem respiration (ER). They are important processes of river ecosystem metabolism (REM) and regulate nutrient and energy flow [53,54]. The difference between GPP and ER indicates the overall carbon budget in river ecosystems, it is therefore referred to as net ecosystem production (NEP) [55]. Use of NEP can be used to determine autotrophic or heterotrophic conditions (e.g., [56,57,58]). A wide range of environmental factors affects both GPP and ER, primarily including light availability [55,59,60,61], water temperature [62,63,64], aquatic species composition [25,65,66] and nutrient supply, notably stressors associated with the USS [35,67]. Overall, considering complex interactions among hydrology, climate, chemistry, and biology, REM has direct implications for water quality and is an excellent functional indicator of ecological response to environmental change [66,68,69,70]. Measurement of REM, therefore, helps to understand the responses of urban river ecosystems to a range of disturbances [53,54,57].

Urban rivers are unfortunately characterized by nutrient enrichment, and there are many studies on elevated concentrations in urban rivers (e.g., [71,72,73]). There are, however, relatively few studies that have integrated ecological processes that incorporate nutrient fluxes and their involvement in river metabolism [66,74]. Further, as a megacity, Beijing’s waterways have experienced multiple nutrient sources and increased nutrient concentrations, resulting in significant changes in river ecological processes [75,76]. In our study, three urban tributaries of the River Beiyun, Beijing, China, provide us with natural conditions with elevated nutrient concentrations to examine river ecological processes. We aim to answer the following questions: (1) How do nutrient dynamics at the sediment-water interface (e.g., concentrations and fluxes of NO₃⁻, NH₄⁺, PO₄³⁻ and DOC) and metabolic processes vary between rivers and between sites on the same river? (2) What are the seasonal characteristics of nutrient dynamics and metabolic processes? (3) Are there relationships among nutrient dynamics and river metabolism? (4) What are the potential drivers of these ecological processes?

2. Materials and Methods

2.1. Study Area

The River Beiyun is an important river in North China with 56% urban land cover, supporting 72% of Beijing’s population [77]. The air temperatures range from −5–26 °C throughout the year and the mean annual precipitation is 440 mm [78]. Three urbanized tributaries were selected in this study (River Liangshui, River Tonghui and River Wenyu) which are WWTP discharge channels for Beijing ([77]; Figure 1). The effluent standard is based on the “Pollutant Discharge Standards for Urban Wastewater Treatment Plants in China (DB11/890-2012)”, including TN at 15 mg/L, TP at 0.3 mg/L, ammonium (NH₄-N) at 1.5 mg/L, biological oxygen demand (BOD) at 6 mg/L, chemical oxygen demand (COD) at 30 mg/L and solids suspended (SS) at 5 mg/L (http://sthjj.beijing.gov.cn, accessed on 29 June 2012). Two study sites were selected on each river, and and their locations were determined by sampling accessibility. These sites were named according to the initials of the river name (i.e., L1, L2, T1, T2, W1, W2), where 1 and 2 represent upstream and downstream sites (relative location towards to the River Beiyun) respectively.

The River Liangshui is 68 km long with a catchment area of 605.7 km² [79]. The distance between the two sites L1 and L2 was 10 km, of which L1 is located at the effluent outlet from the Xiaohongmen WWTP, which serves 2,415,000 people with a mean daily wastewater treatment of 616,000 m³ (https://www.dowater.com, accessed on 24 September 2014). The water depth at L1 is 0.2–0.4 m and the channel width is 15 m. The water depth at L2 is 0.3–0.6 m and the channel width is 30 m.

The River Tonghui is 20 km long with a catchment area of 258.26 km² [80]. The distance between the two sampling sites T1 and T2 is 6 km, of which T1 is located in an urban reservoir formed by the effluent from the Gaobeidian WWTP, which serves 2,400,000 people with a mean daily wastewater treatment of 1,000,000 m³ (https://www.dowater.com, accessed on 7 June 2008) and is used for industrial and agricultural irrigation water. The water depth at T1 is 0.2–0.3 m and the channel width is 100 m. The water depth at T2 is 0.4–0.7 m and the channel width is 10 m. As a reservoir, T1 is bottomed with cement attaching a large amount of aquatic vegetation, a small amount of sediment and smaller phytoplankton; abundant larger phytoplankton is observed in the channel at T2.

The River Wenyu is 48 km long with a catchment area of 4300 km² [81]. The two sites W1 and W2 with a distance of 2 km are located in the downstream area of the River Wenyu, of which W2 is only 1.5 km away from the main river—River Beiyun. Such a close distance from the main river causes the water quality of this reach to have a direct impact on the water environment of the River Beiyun. The River Ba, a tributary of the River Wenyu, flows into the river between the two sites, containing TN at 2 mg/L, TP at 0.4 mg/L, NH4-N at 2 mg/L, BOD at 10 mg/L, COD at 40 mg/L [82], resulting in the difference of water quality between the two sites. The water depth at W1 is 0.3–0.6 m and the channel width is 200 m. The water depth at W2 is 0.4–0.8 m and the channel width is 150 m. Higher water depths were often observed in the rainy season (i.e., summer).

2.2. Sample Collection

River and sediment samples were collected seasonally from September 2020 to August 2021 (autumn from September to October; winter from November to February; spring from March to May; summer from June to August) only at times of stable discharge. In total, 8 sampling events occurred in each season (Table S1). Water pH, water temperature and DO concentration were measured on-site using a multiparameter probe (Hanna HI98194 multiparameter meter). Water samples were taken within 15 cm below the water surface at each site, and a composite 500 mL sample was collected for laboratory testing. The sampling locations for sediments were the same as those for water samples. Ten patches per 10 cm² of the riverbed were randomly selected and a total of 200 mL of sediments (top 2 cm) were collected with a stainless-steel scoop. Water and sediment samples were stored cool separately during transport to the laboratory.

2.3. Concentration Measurement of NO₃⁻, NH₄⁺, PO₄³⁻, DOC

NO₃⁻ and NH₄⁺ concentrations were measured in the field by portable concentration measuring probes (Vernier Go Direct, Beaverton, OR, USA). The PO₄³⁻ concentration was measured by ion chromatography (IC, Dionex Aquion, Thermo Fisher Scientific, Waltham, MA, USA). The DOC concentration was measured by a total organic carbon (TOC) analyzer (TOC-L, SHIMADZU, Kyoto, Japan). A series of concentration gradient standards were prepared to periodically calibrate the instrument. Water samples were passed through 0.45 μm membrane filters before IC and TOC measurements.

2.4. Flux Measurement of NO₃⁻, NH₄⁺, PO₄³⁻, DOC

Flux experiments were run in September and November 2020 and March and June 2021 (Table S2). In the laboratory, to simulate the sediment-water interface, 20 mL of wet surface sediment and 150 mL of overlying water were mixed into a 500 mL chamber, which was scaled up by the nitrate flux experiment by [29]. At 2.5 min after mixing, the particle settlement was visualized. Experiments were run for 3 min and 10 min after mixing to measure the concentrations of NO₃⁻, NH₄⁺, PO₄³⁻, and DOC to represent the process related to both physical disturbance and biological action [29]. We consider that water extraction at t = 3 min would centrally represent physical disturbances, while the experiment where water was extracted at t = 10 min represented a mix of processes. All experiments were run with 3 replicates. A DO logger (PME MiniDot oxygen logger) was also used to continuously measure the DO concentration in the overlying water during two exposures (0–3 min, 3–10 min).

Nutrient fluxes were derived from Equation (1):

$$f=\frac{\left(C_2-C_1\right) V}{A \left(t_2-t_1\right)}$$

(1)

where C₁ and C₂ refer to the nutrient concentration at times t₁ and t₂, respectively; V is the volume of the overlying water (L); A is the surface area of the sediment surface (m₂); t₂ − t₁ is the time (s) between the subsequent (t₂) and previous (t₁) water extraction. Fluxes are expressed as mg/(m²·s). A positive flux indicates movement from sediments to the overlying water and a negative flux indicates movement from the overlying water to the sediment.

During the summer flux experiment, we conducted an additional experiment to further verify the contribution of microbial action to the flux. During the 0–3 min exposure, experimental 3 replicates were run where we added 1 mL of 1 mg/L of zinc chloride with lethal effect on microorganisms [83]. Subsequently, flux was compared to experiments without zinc chloride (Table S3).

2.5. Equilibrium Flux of NO₃⁻, NH₄⁺, PO₄³⁻, DOC

The following nutrient standard solutions were configured: NO₃⁻ standard solution with concentration gradient of 10–50 mg/L with sodium nitrate; NH₄⁺ standard solution with concentration gradient of 0.5–2.0 mg/L with ammonium chloride; PO₄³⁻ standard solution with concentration gradient of 0.5–3.0 mg/L with monopotassium phosphate; DOC standard solution with concentration gradient of 5–20 mg/L with potassium biphthalate (Table S4). In each 500 mL chamber, 20 mL of wet surface sediment was weighed before adding to the chamber, and then 150 mL of standard solution with the above concentration gradient was added. Nutrient concentrations of the overlying solution were measured at t = 3 min. The amount of nutrients released/absorbed by sediments was derived from Equation (2):

$$Q=\frac{(C_0-C_1)V}{W}$$

(2)

where Q is the amount of nutrient released/adsorbed by sediment per unit weight (mg/g); V is the volume of overlying standard solution; W is the weight of wet surface sediment. A positive Q represents the amount of nutrients absorbed by sediment; a negative Q represents the amount of nutrients released from sediment to the overlying standard solution.

Taking the concentration C₀ of the standard solution as the independent variable and Q as the dependent variable, the linear equation Q−C was established for each nutrient (Table S4). When Q = 0, it indicates that the uptake and regeneration flux of nutrients reach equilibrium, and the corresponding nutrient concentration is the equilibrium concentration at the sediment-water interface.

2.6. Sediment Analysis

We collected an additional 150 mL of sediment samples at each sampling site, and loss on ignition [84] was used to quantify sediment organic matter content (OM%). Crucibles were placed in a 120 °C drying oven for 48 h, then their net weights were weighed (m₀). The sediments from each site were divided into 3 equal replicates, spread evenly in a tin foil tray, and dried with crucibles for 48 h. After drying and cooling to room temperature, samples were sieved into coarse sediment and fine sediment using a 2 mm pore size sieve. The coarse and fine sediments were respectively placed in dried crucibles, and the total weight (m₁) of each was measured. m0 was subtracted from m1 to obtain the sediment dry weight (m₂). The crucibles containing the sediment were placed in an oven at 550 °C and burned for 2 h, then placed in a desiccator and cooled to room temperature at which point the weight of each (m₃) was measured (Table S5).

The sediment organic content was derived from Equation (3):

$$OM\%=\frac{m_1-m_3}{m_2}\times 100\%$$

(3)

where m₁ is the total weight before burning in the oven; m₂ is sediment dry weight; m₃ is the total weight after burning in the oven. The grain size of the above fine sediment was measured across all sites using Malvern Mastersizer 2000 Granulometer, which was repeated three times for each site. Samples were classified as either sand (0.063–2 mm), silt (0.004–0.063 mm), or clay (<0.004 mm) [85].

2.7. Measurement of DO and Water Temperature

Under the condition of stable river discharge, lab-calibrated MiniDOT loggers were used to continuously monitor DO concentrations and water temperature across all study sites in September and November 2020 and March and June 2021. DO loggers were placed on the riverbed and measured for 24 h at each site (Table S6). The data collection frequency was set to measure once every 5 min.

2.8. GPP, ER and Reaeration Estimates

We used the open-water, single-station method during stale flows for our REM estimates [55,66,86]. Specifically, the night-time slope modelling (NSM) was applied for GPP, ER and reaeration estimates ([87], Table S7), which is based on the principle of no GPP during nighttime to simplify the oxygen mass balance equation (Equations (4) and (5)) to establish the linear regression between DO change rate and ER. The slope of Equation (4) is the reaeration coefficient k, and the intercept is ER. The nighttime period was selected from the time of sunset to the time when the difference between DO saturation concentration and DO concentration (DO deficit) reached the maximum.

$$\frac{d{C}_{DO}}{dt}=k\left(C_{DOs}-C_{DO}\right)+GPP-ER$$

(4)

$$\frac{d{C}_{DO}}{dt}=k\left(C_{DOs}-C_{DO}\right)-ER$$

(5)

where C_DO is DO concentration (mg/L); C_DOs is DO saturated concentration (mg/L); k is reaeration coefficient (d⁻¹); k(C_DOs-C_DO) is reaeration rate (mg·L⁻¹·d⁻¹).

C_DOs was derived from Equation (6):

$$C_{DOs}=1.4276exp\left[-173.4292+249.6339\left(\frac{100}{T}\right)+143.3483ln\left(\frac{T}{100}\right)-21.8492\left(\frac{T}{100}\right)\right]$$

(6)

where T is mean water temperature (K).

The ER obtained by Equation (5) requires temperature correction (Equation (7)), which is related to Arrhenius activation energy E_a, (i.e., energy required for respiration of aquatic organisms, [88]). The calculated ER was used to obtain GPP in Equation (4).

$${\theta }_R=exp\left(\frac{E_a}{k_B{T}_{0}T}\right)$$

(7)

where θ_R = 0.077, the correction coefficient; E_a is Arrhenius activation energy (eV).

The temperature-corrected ER was derived from Equation (8):

$$E{R}^{*}=ER{\mathsf{\theta }}_R^{T-20°C}$$

(8)

The temperature-corrected GPP was derived from Equation (9):

$$GP{P}^{*}=GPP{\mathsf{\theta }}_P^{T-20°C}$$

(9)

where θ_P = 1.066, the correction coefficient.

The temperature corrected reaeration rate was derived from Equation (10):

$$k\left(C_{DOs}-C_{DO}\right){\mathsf{\theta }}_a^{T-20°C}$$

(10)

where θa = 1.024, the correction coefficient.

2.9. Metabolism Data Standardization

The following equations were applied to convert the values of GPP and ER into the absorption rate or release rate of CO₂ [57]:

$$\frac{gC}{L}=\frac{g{O}_2}{L}\times \frac{12}{32}\times \frac{1}{PQ}$$

(11)

$$\frac{gC}{L}=\frac{g{O}_2}{L}\times \frac{12}{32}\times RQ$$

(12)

where PQ is the photosynthetic quotient (the number of moles of O₂ released per mole of CO₂ incorporated) = 1.2; RQ is the respiratory quotient (the number of moles of CO₂ released per mole of O₂ consumed) = 0.85; the term ‘12/32’ is the ratio of the atomic mass of C to the molecular mass of O₂.

2.10. Structural Equation Modelling

We established pathways for structural equation modelling (SEM) [89] to investigate the effect of environmental and water quality related factors on the REM (Figure 2). Arrows point to the affected factors, boxes represent variables measured in the experiment, and ovals represent latent variables that were not measured. The evaluation indexes of model fitting are provided in Table S8b. The SEM was performed in SPSS 21.0 AMOS.

2.11. Data Analysis

For the study of nutrient dynamics, we sought to evaluate nutrient concentrations and flux differences both among our study sites and between the time exposure experiments. To accomplish this, we used a two-way analysis of variance for each measured water chemistry variable and parameter-specific flux estimates. Regarding the study of river metabolism, we used one-way analysis of variance for each estimated metabolic variable to evaluate the temporospatial differences of REM across study rivers. For the analysis of significant differences between non-normally distributed data, we used non-parametric tests (i.e., Games-Howell). For normally distributed data, a paired t-test was used. Additionally, to explore overall relationships among water chemistry, flux values and metabolic estimates, we examined correlations (i.e., Spearman for non-normally distributed data; Pearson for distributed data) across these parameters. All statistics were performed in SPSS 26.0.

3. Results

3.1. Nutrient Dynamics in Urban River Ecosystems

3.1.1. Water Quality of Study Rivers

Nitrate: During the study period, there were significant differences in NO₃⁻ concentrations among the three rivers (F = 12.737, p < 0.001) and a negative correlation between NO₃⁻ and NH₄⁺ concentrations was observed (r² = −0.557, p < 0.001). Specifically, the concentrations at relative downstream sites were significantly higher than at upstream sites across rivers (p < 0.001, Figure 3a). The maximum and minimum mean NO₃⁻ concentrations appeared at the sites T2 and W1, respectively. Among these rivers, the NO₃⁻ concentration in winter was moderately higher than in autumn (p = 0.310) and then decreased continuously in spring (p < 0.001) and summer (p < 0.001). Overall values ranged from 13.8 to 51.9 mg/L (Figure 3b).

Ammonium: During the study period, there were significant differences in NH₄⁺ concentrations among the three rivers (F = 4.760, p = 0.010). The maximum and minimum mean NH4+ concentrations appeared at the sites L1 and W1, respectively. However, no significant difference in the concentration between the two sites of each river was observed (Figure 3a). A positive correlation between NH₄⁺ concentrations and water temperatures was observed (r² = 0.707, p < 0.001,

Table 1. Seasonal variation of DO concentrations and river metabolic variables across study sites.

Site	Season	DO (mg/L)	Reaeration (mgO2/(L·d))	GPP (mgO2/(L·d))	ER (mgO2/(L·d))	NEP (mgO2/(L·d))	Activation Energy (eV)
L1	Autumn	7.0 (4.3–10.0)	5.3 (−5.3–14.2)	13.3 (0–57.9)	24.7 (21.3–27.9)	−11.4 (−27.9–32.8)	0.57 (0.56–0.57)
	Winter	5.7 (3.9–7.9)	6.2 (3.2–8.8)	4.4 (0–19.2)	8.8 (7.4–11.1)	−4.4 (–9.5–8.8)	0.54 (0.54–0.55)
	Spring	6.7 (5.3–8.7)	13.1 (2.2–19.7)	8.3 (0–28.1)	22.3 (20.7–24.4)	−14.0 (–23.2–5.6)	0.56 (0.56–0.56)
	Summer	8.0 (4.6–11.4)	−1.3 (−6.5–3.9)	20.6 (0–53.2)	16.5 (14.4–19.0)	4.1 (−16.9–35.1)	0.58 (0.57–0.58)
L2	Autumn	7.8 (4.9–11.1)	1.5 (−3.2–5.9)	12.0 (0–33.9)	17.5 (16.3–18.7)	−5.5 (−18.2–15.6)	0.57 (0.56–0.57)
	Winter	6.2 (5.0–7.6)	6.8 (4.6–8.7)	2.3 (0–12.0)	6.9 (6.1–7.8)	−4.6 (−7.1–5.1)	0.54 (0.54–0.55)
	Spring	7.2 (5.5–9.0)	2.2 (−0.7–4.7)	8.3 (0–24.7)	13.6 (12.8–14.5)	−5.3 (−14.5–11.0)	0.57 (0.57–0.58)
	Summer	8.9 (7.3–10.7)	−4.0 (−11.6–2.2)	17.1 (0–46.0)	13.7 (12.6–15.1)	3.4 (−13.9–31.0)	0.58 (0.57–0.58)
T1	Autumn	7.0 (4.1–10.5)	2.3 (−2.2–5.6)	9.6 (0–29.6)	15.5 (14.5–17.2)	−5.8 (−15.7–12.5)	0.56 (0.56–0.57)
	Winter	6.9 (5.5–8.6)	5.5 (2.5–7.7)	2.8 (0–10.3)	7.9 (7.3–8.9)	−5.1 (−8.3–1.4)	0.55 (0.55–0.55)
	Spring	8.5 (6.9–10.3)	0.4 (−4.3–4.5)	7.6 (0–25.8)	11.2 (10.3–12.6)	−3.6 (–11.9–14.0)	0.57 (0.56–0.57)
	Summer	9.1 (7.3–11.3)	−3.3 (−10.4–2.0)	13.4 (0–45.3)	12.0 (10.3–14.6)	1.4 (−13.5–31.7)	0.58 (0.57–0.58)
T2	Autumn	8.4 (6.8–10.3)	1.5 (−4.1–5.6)	8.4 (0–30.5)	11.8 (11.0–13.0)	−3.4 (−12.2–17.6)	0.56 (0.56–0.56)
	Winter	6.8 (5.6–8.4)	4.6 (2.2–6.3)	1.9 (0–9.9)	6.9 (6.4–7.8)	−5.0 (−7.3–2.1)	0.55 (0.55–0.55)
	Spring	8.5 (6.7–10.5)	0.02 (−1.3–1.3)	8.5 (0–23.9)	11.1 (10.7–12.1)	−2.6 (−11.1–12.3)	0.57 (0.57–0.57)
	Summer	10.3 (8.7–12.1)	−8.6 (−16.3–1.6)	19.7 (0–51.0)	8.6 (7.5–10.0)	11.1 (−8.8–41.2)	0.57 (0.57–0.58)
W1	Autumn	7.8 (6.5–9.7)	1.5 (−2.1–3.8)	6.1 (0–25.3)	8.2 (7.6–9.5)	−2.1 (−8.8–15.8)	0.56 (0.56–0.57)
	Winter	6.1 (4.9–7.5)	6.9 (4.6–8.7)	1.9 (0–7.4)	6.4 (5.8–7.5)	−4.5 (−7.0–0.8)	0.55 (0.54–0.55)
	Spring	8.0 (6.6–9.6)	1.7 (−0.4–3.5)	6.0 (0–16.2)	8.0 (7.3–9.0)	−1.9 (−8.5–8.1)	0.56 (0.55–0.56)
	Summer	9.8 (8.2–11.4)	−3.3 (−7.1–0.3)	15.1 (0–42.2)	10.4 (8.7–12.6)	4.7 (−10.7–31.1)	0.57 (0.57–0.58)
W2	Autumn	7.7 (6.2–9.8)	2.0 (−2.6–5.1)	5.7 (0–16.6)	10.6 (9.5–12.9)	−4.9 (−11.4–5.5)	0.56 (0.56–0.57)
	Winter	6.0 (4.8–7.7)	7.8 (5.2–9.7)	2.8 (0–9.8)	7.2 (6.4–8.9)	–4.4 (−7.3–2.7)	0.54 (0.54–0.55)
	Spring	8.1 (6.9–9.4)	2.2 (−0.2–4.3)	5.2 (0–15.2)	7.8 (7.1–8.6)	−2.6 (−8.4–6.6)	0.56 (0.55–0.56)
	Summer	9.3 (8.0–11.0)	5.7 (−14.8–1.4)	15.7 (0–50.7)	7.7 (6.5–9.4)	8.0 (−8.1–41.5)	0.57 (0.57–0.58)

), resulting in the lowest concentration in winter and the highest concentration in summer, with a range of 0.3–3.4 mg/L (Figure 3b). There were significant differences in the NH₄⁺ concentration (p < 0.001) and water temperature (p < 0.001) between seasons.

Phosphate: During the study period, there were significant differences in the PO₄³⁻ concentrations among the three rivers with a range of 0.2–1.5 mg/L (F = 32.916, p < 0.001), and significant differences between the two sites of each river (River Liangshui: p = 0.001; River Tonghui: p = 0.027; River Wenyu: p = 0.038, Figure 3a). The maximum and minimum mean PO₄³⁻ concentrations appeared at W2 and T1 sites, respectively. A negative correlation between PO₄³⁻ and DO concentrations were found (r² = −0.321, p < 0.001), resulting in the highest concentration in winter with the lowest DO. However, the lowest PO₄³⁻ concentrations were recorded in spring (Figure 3b). There were significant differences in the PO₄³⁻ concentrations (p < 0.001 from autumn to winter; p < 0.001 from winter to spring; p = 0.023 from spring to summer) and the DO concentrations (p < 0.001 from autumn to winter; p = 0.001 from winter to spring; p < 0.001 from spring to summer) between seasons.

Dissolved organic carbon: During the study period, there was no difference in the DOC concentration among the three rivers with a range of 4.3–16.7 mg/L, and no significant difference in the concertation between the two sites of each river was observed (Figure 3a). The maximum and minimum mean DOC concentrations appeared at L1 and W1, respectively. Among these rivers, the DOC concentration in winter was significantly higher than in autumn (p = 0.014), and then decreased significantly in spring (p < 0.001) but slightly in summer (Figure 3b). Across all sites, DOC was positively correlated with NO₃⁻ (r² = 0.669, p < 0.001) and PO₄³⁻ (r² = 0.518, p < 0.001) concentrations and negatively correlated with NH₄⁺ concentrations (r² = −0.757, p < 0.001).

3.1.2. Sediment Analysis

Sediments used for the flux experiments varied both among sites. Generally, fine sediments (particle size < 2 mm) contained more organic matter than coarse sediments (particle size > 2 mm). Sediments of all sites were mainly composed of silt, with the highest proportion of fine sediment at L2 except in summer, when the highest proportion of fine sediment appeared at W2 (Table S9). Organic matter (OM%) of sediments at all sites was the highest in summer. OM% at T1 and T2 was the lowest in autumn, while at other sites minimum values were observed in winter.

3.1.3. Nutrient Fluxes at the Sediment-Water Interface

Nitrate and ammonium: Across all sites, NO₃⁻ fluxes in both exposure periods ranged from −0.31 to +0.30 mg/(m²·s) and NH₄⁺ ranged from −0.01 to +0.05 mg/(m²·s) (Figure 4). There were no significant variations in the NO₃⁻ flux between the three rivers, but there were significant differences in flux between the sites on the River Wenyu (p = 0.002). Significant seasonal variations were observed in the River Wenyu (F = 4.129, p = 0.032) but again not for other rivers. No significant spatial differences in the NH₄⁺ flux between the three rivers or between sites on each river were found. NH₄⁺ flux generally presented as regeneration, but no significant seasonal variations were observed among the study rivers. Both NO₃⁻ and NH₄⁺ flux rates were slower in the 3–10 min exposure compared with the 0–3 min exposure.

Phosphate: PO₄³⁻ fluxes in both exposure periods ranged from −0.01 to +0.01 mg/(m²·s) (Figure 5). For the River Liangshui, L2 showed the maximum regeneration flux in autumn and L1 showed the maximum uptake flux in winter, but no overall significant seasonal variations were recorded. The two sites on the River Tonghui only presented regeneration fluxes in the 0–3 min exposure in autumn. Similarly, seasonal variations did not vary significantly. The two sites on the River Wenyu presented uptake fluxes in the 0–3 min exposure in all seasons, but again no significant seasonal variations were observed. Spatially, no significant differences among rivers or between sites on each river were observed.

Dissolved organic carbon: DOC fluxes in both exposure periods ranged from −0.04 to +0.13 mg/(m²·s) (Figure 5). In the 0–3 min exposure, DOC at the sites on the River Liangshui and the River Tonghui was released from sediments except for L1 in autumn. In the 3–10 min exposure, DOC at the sites on the River Wenyu showed uptake flux except for W1 in autumn. The maximum uptake and regeneration fluxes occurred at W2 in winter and at T1 in summer, respectively. There were no overall significant variations in the DOC flux across the three rivers or between the sites on each river. Further, no significant seasonal variations were observed in each river.

DO flux: In both exposures, DO was consumed across all sites with the consumption flux ranging from 0.01 to 0.05 mg/(m²·s) (Figure S1). The DO consumption flux was positively correlated with regeneration fluxes of NH₄⁺ (p < 0.001) and DOC (p < 0.001), but no correlation was observed with NO₃⁻ or PO₄³⁻ flux; only the initial PO₄³⁻ concentration was positively correlated with PO₄³⁻ regeneration flux (p = 0.002), indicating that in the overlying water DO and initial nutrient concentrations did not independently affect nutrient fluxes.

When the water-specific nutrient concentration increased/decreased seasonally, an increase in the equilibrium concentration strongly indicated an increase in the sediment-specific nutrient concentration, while a decrease in the equilibrium concentration indicated a decrease in the sediment-specific nutrient concentration (Table S10). When the water-specific nutrient concentration was higher than the sediment-specific concentration, the nutrient showed uptake flux, otherwise, a regeneration flux was observed.

In the summer flux experiment, only a few sites with decreased nutrient fluxes and DO consumption fluxes were observed in the experimental group containing 1 mg/L zinc chloride (Figure S2). There was no significant difference between the experimental group and the control group, indicating that the main factor affecting nutrient flux in the 0–3 min exposure was physical disturbance. This result supports the contention regarding the roles of physical process and biological process during our two experimental exposure periods, respectively.

3.2. Metabolic Patterns of Urban River Ecosystems

3.2.1. DO Profiles among Study Sites

The minimum and maximum mean DO concentrations occurred at all sites in winter and summer, respectively (Figure S3). Significant seasonal differences in DO were observed among the study sites (p < 0.001) except T1, where there was no significant difference between autumn and winter (Table 1). Regarding diurnal DO, the maximum DO concentration occurred between 1 and 5 pm during the four seasons and decreased gradually at night. The DO concentration increased again when the maximum DO deficit was reached. The maximum DO deficit appeared earliest in winter, followed by spring and autumn, and latest in summer.

3.2.2. GPP, ER and NEP Estimates

The minimum mean GPP appeared in winter and decreased from noon; the maximum mean GPP appeared in summer and peaked between 1 and 3 pm (Figure 6). The mean GPP in autumn was higher than in spring except for T2, while the maximum GPP at T2 in autumn was higher than spring values. Compared with autumn, GPP at L2 and T2 decreased significantly in winter (L2: p = 0.011; T2: p = 0.028) then increased significantly in spring (L2: p = 0.048; T2: p = 0.020); GPP at W2 increased significantly from spring to summer (p = 0.044). However, no significant differences were observed in GPP for all rivers spatially.

ER among the three rivers was the lowest in winter and the highest in autumn except for W1, where the highest ER was observed in summer (Figure 7). Significant seasonal variations at all sites were recorded from autumn to spring (p < 0.001, except p = 0.006 for W2 from winter to spring); in summer, ER at L1, T2 and W1 still increased significantly (p < 0.001). Different from GPP, significant differences in ER were observed between the two sites of each river (p < 0.001, except River Tonghui and River Wenyu in spring). Further, significant spatial differences in ER were also observed across the three rivers except River Liangshui and River Tonghui in winter (p = 0.002 between River Tonghui and River Wenyu in summer; p < 0.001 among the three rivers in other seasons). There was a positive correlation between GPP and ER across the three rivers during the study period (r² = 0.462, p < 0.001).

In each season, NEP decreased to the lowest value within one hour after sunset and then gradually increased (Figure 8). The maximum mean NEP across all sites appeared in summer and reached the peak between 10 am–2 pm. The minimum NEP at L1 was observed in spring, while the minimum values of other sites appeared in autumn when diurnal NEP peaked between 10 am–12 pm. The summertime positive mean NEP indicated that these river ecosystems were overall sinks of atmospheric CO₂. For these periods mean absorption rates ranged from 1.6 to 12.7 gCO₂/(m³·d). At other times, the mean NEP values were negative indicating these sites were sources of CO₂ with the mean release rate ranging from 2.2 to 16.0 gCO₂/(m³·d). The springtime NEP in the River Tonghui was significantly different from the River Liangshui (p = 0.006) and the River Wenyu (p < 0.001). Also, there was a significant difference between the two sites of the River Liangshui in spring (p = 0.005).

3.2.3. Reaeration and Activation Energy Estimates

In each season, the maximum reaeration rate across all sites appeared within one hour before sunrise and gradually decreased with DO accumulation in the GPP process during the daytime (Figure S4). The seasonal reaeration rates were the lowest in summer among all sites, especially for T2, the negative values indicated river DO saturation. Except L1, the highest reaeration rates were observed in winter. For L1, the lowest NEP in spring contributed to the highest reaeration rate in this season.

Significant seasonal variations in the reaeration rate were observed in all rivers between each season (p = 0.031 for L1 from autumn to winter; p = 0.006 for T1 from spring to summer; p = 0.001 for T2 from autumn to winter; p < 0.001 for other sites in other seasons). In winter and spring, there were significant spatial differences among each river (p < 0.001, except p = 0.007 between River Liangshui and River Wenyu in winter), with the highest mean rate appearing at W2 in winter and L1 in spring, which were 7.8 mgO₂/(L·d) and 13.1 mgO₂/(L·d), respectively. Across all rivers, significant differences were also observed between each site in summer (p < 0.001).

Among all sites, the activation energy (E_a) required for ecosystem respiration was in the range of 0.54–0.58 eV during the study period. The minimum and maximum mean E_a of each river were recorded in winter and summer respectively and significant temporal differences were observed for all sites and among all seasons (p < 0.001). In terms of the spatial distribution of E_a, similar to ER, significant differences were observed among the three rivers in autumn, spring and summer (p = 0.019 between River Tonghui and River Wenyu in summer; p < 0.001 between other rivers in the three seasons). However, the significant difference in E_a at the two sites on each river was not always observed, especially for the River Wenyu, there was no significant difference in all seasons.

3.2.4. SEM Pathway Analysis

Structural equation model (SEM) revealed the selected effects of environmental and water quality related factors on river ecosystem metabolism. Increased concentrations of NH₄⁺ and PO₄³⁻ significantly affected bioactivity then increased GPP and ER (Figure 9). Unsurprisingly, water temperature was significantly related to bioactivity which was positively correlated with GPP (r² = 0.89, p < 0.001) and ER (r² = 0.58, p = 0.003). Elevated water temperature which reduces the solubility of oxygen was observed as a negative correlation with DO (r² = −0.90, p = 0.020). ER was negatively correlated with DO concentrations, however, there was not a significant correlation between GPP and DO, suggesting that ER played a key role in consuming DO. Seasonal significant variations of DO determined the reaeration process, resulting in significant seasonal differences in the reaeration rates across all sites being recorded.

4. Discussion

4.1. Nutrient Dynamics across the Sediment-Water Interface

Our results indicated that discharged WWTP effluents and tributaries rich in nutrients resulted in elevated nutrient concentrations and reduced DO concentrations at several sampling sites (i.e., L1, T1, W2; Figure 3). We found that the sediments at T1 and L1 were the largest source and sink of NH₄⁺ in summer, respectively. Further, sites L1 and W2 were the largest source and sink of PO₄³⁻ in autumn, respectively (Figure 4 and Figure 5). Specifically, we found that the rate and direction of nutrient flux at the sediment-water interface were affected by nutrient equilibrium concentration, DO conditions in the overlying water and likely patch-scale microbial activity in the sediment.

The concentration changes across the sediment-water interface are influenced by a combination of physical processes and biogeochemical processes. In summer when we suspect there would be the highest microbial activity, no significant differences were observed in the 0–3 min exposure for either nutrient fluxes or DO consumption between the control group and the experimental group with ZnCl₂. These observations support that in the first exposure period, physical disturbances caused by water movement and sedimentation were the main factors affecting nutrient fluxes. For the subsequent experimental period (i.e., 3–10 min exposure), microbial processes were the likely influencing factors but further studies are needed to quantify these factors.

Generally, the role of sediments as NO₃⁻ producers/consumers was accompanied by the role as NH₄⁺ consumers/producers (Figure 4), which was related to nitrogen transformations (e.g., nitrification and dissimilatory nitrate reduction to ammonium (DNRA)). The transformations across the sediment-water interface were attributed to the seasonal variations of NO₃⁻ and NH₄⁺ concentrations, water temperature and DO [90,91]. For PO₄³⁻, insoluble or soluble matter in the sediment affected by DO and water pH interacted with minerals and associated adsorption sites. Therefore variation in the adsorption capacity of sediments was likely responsible for observed differences [38,40] (Figure 5). DOC dynamics related to its roles as a biological energy source and respiration substrate were likely affected by major processes related to water column ecosystem metabolism and the overall content of organic matter in sediments [24,25] (Figure 5; Table S9).

4.1.1. Mechanisms Driving NO₃⁻ and NH₄⁺ Concentrations and Fluxes

Among the study rivers, NO₃⁻ concentrations were highest in winter and the lowest in summer (Figure 3b). On the contrary, NH₄⁺ concentrations were highest in summer and lowest in winter. After water and sediment samples were mixed, disturbances across the sediment-water interface for the 0–3 min exposure caused NO₃⁻ and NH₄⁺ to migrate along the concentration gradients. In winter, the lowest NO₃⁻ equilibrium concentration (26.4 mg/L) and the highest water concentration (43.3 mg/L) were found at site L2 and resulted in the maximum NO₃⁻ uptake flux rate (0.31 mg/(m²·s)) observed in the 0–3 min exposure (Figure 4). For site L1 (i.e., highest NO₃⁻ equilibrium concentration among all sites in winter (45.5 mg/L) but lower water concentration than L2), we observed NO₃⁻ regeneration. For NH₄⁺ in winter, the water concentration at L2 was lower than at L1, while the equilibrium concentration at L2 was the highest among all sites (3.9 mg/L). This likely indicated higher sediment NH₄⁺ concentrations here which resulted in a greater concentration gradient and an overall higher regeneration rate (0.04 mg/(m²·s); Figure 4).

In summer the equilibrium NO₃⁻ concentration at L1 dropped to the lowest in the four seasons (20.9 mg/L) and was lower than the mean water concentration (26.2 mg/L). These conditions indicated lower sediment concentrations than in the overlying water and contributed to NO₃⁻ uptake in the 0–3 min exposure. The NO₃⁻ concentration change at L2 was similar to that at L1, hence the flux still presented uptake in this exposure. The likely decreased concentration gradient also led to a decreased flux rate from 0.09 mg/(m²·s) in spring to 0.05 mg/(m²·s) in summer. For NH₄⁺ in summer, water concentrations at L1 were lower than at L2, while the equilibrium concentration at L1 was higher. This suggests that higher sediment NH₄⁺ concentrations at L1 resulted in a greater concentration gradient of NH₄⁺ and this was observed as regeneration (0.04 mg/(m²·s)). For spring and autumn, as well as other study rivers, seasonal variations in the water concentrations and equilibrium concentrations of NO₃⁻ and NH₄⁺ caused flux transformation to vary between uptake and regeneration (Figure 4). However, the highest water concentration or equilibrium concentration did not correspond to the maximum flux rates.

In terms of biological processes (i.e., 3–10 min exposure), the movement of water transferred NO₃⁻ into the sediment and likely made it bioavailable to microbes associated with sediment surfaces. The DO consumption during both exposures (Figure S1) favored the occurrence of the dissimilatory nitrate reduction to ammonium (DNRA), so that NO₃⁻ was transformed into NH₄⁺ and returned to the overlying water along these concentration gradients [92]. This was consistent with both NO₃⁻ uptake and NH₄⁺ regeneration fluxes shown at some sites (e.g., L1 and L2 in winter; T1 and T2 in summer; Figure 4). However, when the overlying water could not provide sufficient NO₃⁻ DNRA were likely inhibited resulting in decreased NH₄⁺ regeneration fluxes and overall sediment-associated uptake.

Seasonal variations of NO₃⁻ and NH₄⁺ fluxes were determined by a combination of physical and biological factors. There was no significant seasonal difference in the NO₃⁻ equilibrium concentration among the three rivers which resulted in a lack of seasonal difference in the NO₃⁻ flux in the 0–3 min exposure. On the contrary, significant seasonal differences in the NH₄⁺ equilibrium concentration (F = 30.272, p < 0.001) led to the observed seasonal difference in NH₄⁺ flux (F = 34.857, p < 0.001). For the 3–10 min exposure, seasonal differences in water temperature and DO (p < 0.001) affected NH₄⁺ fluxes. We assume this was due to microbial activity affecting the N-cycle and in turn, these seasonal differences affected bioactivity (F = 17.276, p < 0.001). Regarding NO₃⁻, significant seasonal differences in flux associated with bioactivity were also observed (F = 8.723, p = 0.031). We suspect this significant variation reflects NH₄⁺ flux participating in overall DNRA reactions.

Concentrations of NO₃⁻ and NH₄⁺ were significantly different among the study rivers (F = 12.737, p < 0.001; F = 4.760, p = 0.010, respectively). The NH₄⁺ concentrations at the sites affected by the wastewater effluents (i.e., L1, T1) and the site affected by the tributary inflow (i.e., W2) were higher than at another site on each river. On the one hand, the effluents with elevated water temperature and organic matter promoted the mineralization of organic nitrogen by microbes likely promoting higher NH₄⁺ concentrations [92]. On the other hand, high concentrations of NH₄⁺ (>2 mg/L) of the River Ba flowing into the River Wenyu yielded higher NH₄⁺ concentrations at W2. For NO₃⁻, concentrations were negatively correlated with NH₄⁺ concentrations (r² = −0.557, p < 0.001) and NO₃⁻ at L1 and T1 were lower than at L2 and T2, respectively. [93] pointed out that NO₃⁻ tend to increase downstream with increasing river discharge due to reduced water residence time which could explain our observations. Mean NO₃⁻ concentrations across all rivers were lowest in summer. We suspect that this observation is linked to NO₃⁻ uptake supporting biological activity in summer (e.g., phytoplankton). This idea is further supported by the fact that the highest mean NO₃⁻ concentrations were observed in winter. Mean NO₃⁻ concentration in autumn was higher than in spring, and the increased seasonal precipitation in North China is a likely factor because of the role of surface runoff as a non-point source of pollution [14]. Although increased precipitation and runoff can also dilute NO₃⁻ concentrations, pollution associated with surface runoff can be significantly stronger than dilution [94]. However, this needs further investigation for the study basins.

4.1.2. Mechanisms Driving Po₄³⁻ Concentrations and Fluxes

In our study, higher PO₄³⁻ concentrations at L1 and T1 were attributable to WWTP effluents, while the higher concentration at W2 was related to the tributary inflow (Figure 3). However, contributions from PO₄³⁻ released from sediments were less which was reflected in uptake flux and low regeneration rates (Figure 5). In addition to the PO₄³⁻ flux determined by the equilibrium concentration, abundant biofilms at T1 and W2 likely reduced the interaction between PO₄³⁻ and sediment adsorption sites, thereby reducing regeneration from sediments [52,95].

We observed that PO₄³⁻ concentrations were negatively correlated with DO (r² = −0.321, p < 0.001). This finding was associated with the highest mean PO₄³⁻ concentration in winter when the lowest mean DO concentration occurred (Figure 3b). Also, the reduced PO₄³⁻ absorption by weaker bioactivity in winter contributed to increased PO₄³⁻ concentrations [10,46,96]. DO conditions affect electron receptors in sediments and under reducing conditions the formation of soluble substances (e.g., Fe(OH)₂) reduces the adsorption sites resulting in increased water-specific PO₄³⁻ concentrations [39,41]. However, the highest mean DO concentration in autumn did not lead to the lowest PO₄³⁻ concentration and this could be due to increased seasonal precipitation carrying PO₄³⁻ from urban impervious surfaces into rivers [97].

There was a positive correlation between PO₄³⁻ and water pH in the River Tonghui (r² = 0.257, p = 0.040). These study reaches were alkaline (pH = 7.5–8.9), causing OH⁻ and PO₄³⁻ to compete for adsorption sites on the sediment surface and exchange with PO₄³⁻, consequently elevating PO₄³⁻ concentration in the water [98]. However, for the other two study rivers, positive correlations were not observed as sediment particle size seemed to also affect the P-dynamics. For example, both the sediments at L2 and W2 contained a high proportion of fine particles (Table S9) which were easy to resuspend and could cause PO₄³⁻ to re-enter the overlying water [99].

In the 0–3 min exposure, PO₄³⁻ concentrations in water across the study sites were both higher/lower than the equilibrium concentration resulting in both uptake/regeneration fluxes. In the 3–10 min exposure, decreased DO facilitated releases of PO₄³⁻ from sediments and regeneration fluxes or decreased uptake fluxes were observed. Significant seasonal variations in DO and PO₄³⁻ concentrations in water were recorded, but again there were no significant seasonal differences in the PO₄³⁻ flux, which may be attributed to the small PO₄³⁻ flux rate found at lower concentrations when compared with results from the other nutrients.

4.1.3. Mechanisms Driving DOC Concentrations and Fluxes

DOC concentrations increased across all sites in autumn (Figure 3b) and were likely linked to urban impervious surfaces and stormwater runoff [100]. Additionally, tree leaves and other terrestrial organic matter are other DOC sources [28]. However, the highest mean concentration of DOC was found in winter when decreased bioactivity reduced biological uptake [24,25].

The DOC equilibrium concentration determined the DOC flux in the 0–3 min exposure. The content of organic matter in sediments varied little with seasons (Table S9) but there were significant differences in the DOC equilibrium concentration between winter and spring across the three rivers (p = 0.003). These differences were related to the significant difference in DOC flux between these two seasons in both exposure times (p = 0.045). DOC uptake fluxes or decreased regeneration fluxes were observed at all sites in the 3–10 min exposure and this was likely because microorganisms increased the saturation index of biogenic carbonate minerals through their metabolic activities [101,102]. The fine sediments across all sites were mainly composed of silt (Table S9). These sediments have a high surface area available for microbial adhesion and colonization resulting in high microbial biomass and diversity [93]. Therefore, we suggest that all sites may have similar microbial biomass and diversity and this may account for the similar DOC flux among all sites. More detailed studies to confirm these processes are needed.

4.1.4. Flux Predicted by Nutrient Concentrations and Water Quality

We used forward selection stepwise multiple regression to evaluate if measured flux rates could be predicted from ambient river water quality (

Table 2. Stepwise regression (forward selection) for each nutrient fluxes with nutrient concentrations and water temperature.

	NO3− Concentration	NH4+ Concentration	PO43− Concentration	DOC Concentration	Water Temperature	DO
NO3− flux F = 4.368 p = 0.042		𝛽 = −0.525 R2 = 0.292 P = 0.042
NH4+ flux F = 28.185 p < 0.001					𝛽 = 0.612 r2 = 0.612 p < 0.001
PO43− flux F = 6.194 P = 0.016			𝛽 = −0.622 r2 = 0.341 p = 0.016			𝛽 = −0.662 r2 = 0.453 p = 0.001
DOC flux F = 12.158 p = 0.001				𝛽 = −0.277 r2 = 0.113 p = 0.043

Note: Data involved in the regression analysis contains both exposure periods. The final model summary statistics are reported for each flux (below label). Individual parameter fit summary statistics for selected models are given in grid cells. The blacked-out grid cells indicate that these parameters were not selected by the final model. 𝛽 is the standardized regression coefficient; r² is the proportion of variance for a dependent variable that is explained by an independent variable; p is the significance value for each parameter reported.

). Only NH₄⁺ concentration explained NO₃⁻ flux differences suggesting that nitrification plays a strong role within the N-cycle which is similar to results from other rivers with high nutrient concentrations (e.g., [17,103]). For NH₄⁺ flux, only water temperature explained overall differences, pointing towards microbial activity associated with assimilative nitrogen uptake as the driving factor [104]. For PO₄³⁻ flux, PO₄³⁻ and DO concentrations explained overall differences. This suggests that DO conditions affected electron receptors in sediments, thereby affecting the PO₄³⁻ adsorption sites. However, this negative relationship suggests that increased concentrations in river water may saturate uptake [39,41]. Finally, DOC flux was explained by DOC concentration. The negative beta coefficient again suggests that increased DOC leads to regeneration rather than uptake of DOC, demonstrated in the shift of DOC flux from uptake to regeneration. Evaluating each of these models provides an indication of the relative importance of physical and biological processes which drive the dynamics of nutrients at the sediment-water interface in urban rivers.

4.2. Metabolic Patterns of Urban River Ecosystems

River ecosystem metabolism (REM) is well recognized as a holistic measure of ecological function given its integration of complex biological and physical interactions, which integrates gross primary production (GPP) and ecosystem respiration (ER) [54,67]. Pollution stressors (e.g., excessive nutrients) and water temperature drive variations in REM patterns [105]. Urban rivers often receive WWTP effluents and tributary inflows which can alter nutrient concentrations and modify temperature regimes of receiving water bodies, thus affecting REM [18,63,68]. Our SEM modelling revealed that NH₄⁺ and PO₄³⁻ were key nutrients stimulating aquatic bioactivity, thus affecting GPP and ER, and ER significantly consumed DO. In our study, GPP was similar at each site, but ER was significantly higher at the sites affected by the WWTP effluents or the tributary inflow. In summer, the mean NEP at each site was greater than 0, indicating that each study reach was autotrophic and overall carbon sinks with a mean absorption rate ranging from 1.6 to 12.7 gCO₂/(m³·d). However, in other seasons negative NEP values indicating a heterotrophic state of these river reaches and overall sources of CO₂ with a mean release rate of 2.2–16.0 gCO₂/(m³·d) across study sites. This work highlights the major environmental factors driving REM and allows for a holistic understanding of the variation of these processes which suggest options for in-stream management.

4.2.1. DO Profiles across Study Sites

The continuous monitoring of diurnal DO concentrations provided rich information on DO conditions such as the mean and range of DO, especially low DO, which is a useful indicator of habitat suitability for aquatic organisms [106,107,108]. We found that the magnitude of diurnal DO in summer and autumn was higher than in spring and winter due to high rates of photosynthesis and organic matter respiration (Table 1). These conditions can cause physiological stress in aquatic organisms [106].

Regarding spatial distributions of DO, the mean DO concentration at L1, T1 and W2 affected by the WWTP effluents or the tributary inflow was lower than at another site of the river, and the DOC concentration at these affected sites was higher, indicating that organic matter decomposition and respiration contributed to DO consumption (Figure 3; e.g., the maximum DO deficit at L1 in the four seasons).

4.2.2. Factors Driving GPP and ER

In this study, the WWTP effluents and the tributary inflow were used as natural experiments to investigate the effects of elevated nutrient concentrations and water temperature on urban river metabolism. We found similar GPP across our study sites but higher ER in the effluent/tributary-affected reaches (Table 1). [109] proposed criteria to assess the functioning impact on river ecosystems according to GPP and ER values (

Table 3. Mean and seasonal range of GPP and ER (both were multiplied by water depth) across the study sites during the study period. The minimum and maximum GPP appeared in winter and summer, respectively. The minimum and maximum ER appeared in winter and autumn (or summer at W1), respectively. Impacts relate to the degree of impaired ecosystem functioning.

Site	GPP (Min Mean–Max Mean) (gO2·m−2·d−1)	Impact	ER (Min Mean–Max Mean) (gO2·m−2·d−1)	Impact
L1	0.9–8.2	No to severe effect	1.8–7.4	Mild to severe effect
L2	0.7–10.3	No to severe effect	2.1–8.8	Mild to severe effect
T1	0.6–4.0	Mild to severe effect	1.6–4.6	Severe effect
T2	0.8–13.8	No to severe effect	2.8–7.1	Mild to severe effect
W1	0.6–9.1	No to severe effect	1.9–6.2	Mild to severe effect
W2	1.1–12.5	No to severe effect	2.9–7.4	Mild to severe effect

). Based on this, the ecological function of the WWTP effluent reservoir where T1 is located was mildly to severely impaired; the ecological function at other sites changed from no functional effect in summer or autumn to severe impairment in winter. Although GPP and ER were measured at specific locations, they can represent the metabolic patterns of the studied reaches and can be used to assess ecological functions, because of the natural movement and mixing of water in a river [109]. Overall, this suggests that both the WWTP effluents and the tributary rich in nutrients have negative impacts on river ecological function and health for urbanized catchments.

Our SEM model revealed that NH₄⁺ and PO₄³⁻ were key nutrients stimulating aquatic bioactivity (Figure 9), which was directly associated with the promotion of growth and reproduction of aquatic organisms (e.g., macrophytes) and increased both GPP and ER [1,2]. For example, the mean GPP and ER at L1 with effluent outfall were higher than at L2 in the four seasons (Table 1). [74] referred to the increase in GPP and ER in heterotrophic rivers as the “subsidy” effect of WWTP effluents. As a result of this, a positive correlation between GPP and ER was found (r² = 0.262, p < 0.001), which has also been widely reported by many studies (e.g., [110,111,112]). However, fine particulate organic matter (FPOM) and DOC associated with WWTP effluents or tributary inflows can be an important energy input for supporting in-situ ER [74,113]. Additionally, GPP is directly controlled by light availability [66,114]. Therefore, the increase for ER and GPP in reaches with elevated nutrients is not always consistent [55], which explains the inconsistency between the increased significance of GPP and ER across the study sites. Furthermore, the change of light availability led to the diurnal and seasonal variations of GPP as we observed that the mean GPP was the highest in summer across all sites with the peak between 1–3 pm, followed by autumn and spring and was the lowest in winter with the peak at noon (Figure 6).

Regarding the River Wenyu, GPP at W2 with higher nutrient concentrations was not always higher than at W1 (i.e., in autumn and spring) because in addition to nutrient supply, water temperature, turbidity and light availability also likely led to GPP differences [55,115]. We suspect that this is especially the case for light availability as the water at W2 was deeper than at W1. This difference in depth resulted in more attenuation of light at the riverbed depth at W2 [116]. These conditions reduced photosynthesis at the riverbed and caused lower observed DO concentrations. There was a similar situation in the River Tonghui, where at T1 (e.g., a site where a pool is formed by the WWTP effluent) GPP was not always higher than at T2 (i.e., in spring and summer). A potential explanation for this observation is related to algal senescence at T1 [117].

ER was positively correlated with water temperature (r² = 0.567, p < 0.001) due to temperature-related enzyme activity involved in respiration [63,64]. The SEM model demonstrated that bioactivity was significantly affected by water temperature. As a result, due to the relatively warm WWTP effluents, ER at L1 and T1 was significantly higher than at L2 and T2, respectively (p < 0.001). Also, the lowest mean ER was observed in winter when the water temperature was the lowest. The highest mean ER occurred in autumn (Figure 7), and we suggest that this was associated with higher DOC concentrations providing additional energy sources serving as respiratory substrates [24,25]. Of course, ER cannot increase indefinitely with increasing temperature due to factors like substrate limitation and/or protein degradation [25,118]. As a result, with the continuous DO consumption by ER and the potentially limited carbon supply due to the metabolism of carbon compounds, the maximum DO deficit was reached before sunrise and the initiation of photosynthesis [87]. In some systems, this can produce anoxic or hypoxic periods [119], but this was not the case during our study.

Regarding activation energy (E_a) associated with ER, the E_a at each site changed significantly in each season (Table 1). Higher E_a suggests that ER is more sensitive to temperatures and more vulnerable to climate change [120]. Ref. [118] pointed out that catabolic processes are more sensitive to changes in temperature than anabolic processes, and thus, warming should have stronger effects on ER than GPP, suggesting that global warming can reduce stream NEP [121].

The balance of GPP and ER plays an important role in carbon cycling in urban river ecosystems [122]. As we found, high nutrient concentrations significantly increased ER which has implications for CO₂ emissions. Our study reaches released CO₂ at a rate of 2.2–16.0 g/(m³·d) during the study period except for summer, when each reach was autotrophic as overall sinks of atmospheric CO₂ with a mean absorption rate ranging from 1.6 to 12.7 g/(m³·d). Strategies to reduce nutrient concentrations in urban rivers via improved river water quality are paramount for reducing overall ER and thereby lowering the associated atmospheric CO₂ emissions.

4.2.3. Factors Driving Reaeration

At night the lack of oxygen produced by photosynthesis in rivers significantly promoted reaeration, additionally lower water temperatures at night increased the solubility of oxygen [60]. As a result, diurnal reaeration rates were always positive at night and sometimes negative in the daytime (Figure S4) representing input and output of oxygen to the river, respectively. In summer, the highest water temperature and GPP supported high DO concentrations resulting in negative mean reaeration rates at all sites. In other seasons, reaeration supplemented oxygen into the water (Figure 9). Between each season, a significant difference in the reaeration rate at each site was observed which was also attributable to significant seasonal variation in ER. Our SEM model indicated that ER played a decisive role in the DO consumption, causing the significant seasonal difference of DO at each site (Table 1) and thereby enhancing the reaeration rate. For future river restoration, increasing oxygen by promoting reaeration, perhaps via turbulent mixing associated with geomorphic alteration is a viable option [123].

4.3. Research Limitations

The flux chamber experiment is a simple, cheap and easy method and allows us to separate the factors affecting nutrient fluxes into physical and biological processes. However, as chamber experiments, our work cannot capture in situ hydraulic conditions. Although we collected water at different depths, there was still a lack of water advection and hyporheic flow in the chambers, which clearly could alter observed fluxes versus the actual in-stream fluxes. Nevertheless, spatial and seasonal variation of nutrient fluxes across our study still experimentally demonstrated nutrient source and sink dynamics at the sediment-water interface. These data are useful for providing support in considering how to approach environmental protection against sediment degradation.

For river metabolism study, the use of reliable, measuring diurnal and seasonal patterns of river metabolism based on DO change rate from low-cost and reliable data loggers allowed for the estimation of GPP, ER and reaeration rates. This method is an overall estimation of ecosystem metabolism in the study reach; however, it cannot indicate the proportion of GPP and ER in different communities (e.g., phytoplankton and photosynthetic bacteria), which needs further research. Additionally, greater spatial coverage of measurements would help to elucidate reach-specific conditions including overall ecosystem metabolism.

5. Conclusions

This research evaluated urban river ecosystem function by studying flux dynamics of nutrients and river metabolic patterns. Three urban tributaries of the River Beiyun in Beijing each characterized by WWTP effluents and tributary inflows provided conditions with elevated nutrient concentrations to examine river ecological response in terms of nutrient and energy flow. While conditions were significantly different across the study sites, we did find that nutrient concentrations modified biogeochemical flux. We also found that nutrient fluxes were driven by physical processes (i.e., 0–3 min exposure experiments). Due to physical disturbance, nutrients moved along concentration gradients at the sediment-water interface. However, as nutrients were brought into contact with sediment particles, they were absorbed and utilized by microorganisms which speed biogeochemical reactions. Additionally, NO₃⁻ and NH₄⁺ fluxes were seasonally affected by water temperature, DO and presumed bioactivity. We suggest that, especially in summer, this leads to an increased risk of eutrophication. PO₄³⁻ flux was primarily related to the availability of suitable adsorption sites on the sediment surface which related to DO, water pH and organic matter content. Particle size and texture related to sediment resuspension contributed to PO₄³⁻ regeneration and perhaps overall bioactivity (i.e., observed as variations in DOC flux). Regarding PO₄³⁻ and DOC, elevated concentrations lead to regeneration and increase risks of river eutrophication.

Among these nutrients, we suggest that increased NH₄⁺ and PO₄³⁻ significantly stimulated bioactivity promoting GPP and ER which was demonstrated by positive correlations in our study. Our SEM model revealed that ER was related to significant seasonal variations of reaeration rate and DO consumption. Among these rivers, the highest mean GPP was recorded in summer, but the highest mean ER occurred in autumn which is likely linked to increased DOC acting as a carbon source for the heterotrophic processes. Overall higher GPP than ER in summer (NEP > 0) provides evidence to show that the studied river ecosystems were autotrophic; however, in the other seasons, sites were heterotrophic (NEP < 0) and became sources of CO₂.