Algal proliferation risk assessment using Vine Copula-based coupling methods in the South-to-North Water Diversion Project of China

Nong, Xizhi; Zeng, Jun; Ma, Yonghang; Chen, Lihua; Zhang, Chi; Behzadian, Kourosh; Campos, Luiza C.

doi:10.3389/fevo.2023.1193163

ORIGINAL RESEARCH article

Front. Ecol. Evol., 13 June 2023

Sec. Interdisciplinary Climate Studies

Volume 11 - 2023 | https://doi.org/10.3389/fevo.2023.1193163

This article is part of the Research Topic Hydrological and Hydrochemical Processes in Fragile Ecosystems View all 9 articles

Algal proliferation risk assessment using Vine Copula-based coupling methods in the South-to-North Water Diversion Project of China

Xizhi Nong^1,2 $Jun Zeng\n$ Jun Zeng¹ $Yonghang Ma\n$ Yonghang Ma¹ $Lihua Chen\n$ Lihua Chen¹ $Chi Zhang\n$ Chi Zhang²

Kourosh Behzadian^3,4

Luiza C. Campos⁴^*

¹College of Civil Engineering and Architecture, Guangxi University, Nanning, China
²State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan, China
³School of Computing and Engineering, University of West London, London, United Kingdom
⁴Department of Civil, Environmental and Geomatic Engineering, University College London, London, United Kingdom

The Middle Route of the South-to-North Water Diversion Project of China (MRSNWDPC), i.e., the longest inter-basin water diversion project (1,432 km) in the world, has delivered more than 60 billion m³ of water resources to North China and benefiting more than 100 million people since December 2014. However, the abnormal algal proliferation in the main canal under low nutrient background has seriously threatened the water quality safety of this mega project. In this research, 3 years of monitoring data matrix, including water temperature (WT), flow discharge (Q), flow velocity (V), dissolved oxygen (DO), and the algal cell density (ACD), from the main canal of the MRSNWDPC were analyzed. The nonlinear relationships were determined based on multiple regression models, and a composite risk analysis model was constructed by Latin hypercube sampling (LHS) method coupled with Vine Copula function. The impacts of different hydrological and environmental factors on algal proliferation were comprehensively analyzed by Bayesian theory. The results showed that the WT gradually decreased from upstream to downstream, with a narrow range of 16.6–17.4°C, and the annual average concentrations of DO showed a gradual increase from upstream to downstream. The flow velocity of MRSNWDPC had a tendency to increase year by year, and the maximum flow velocity exceeds 0.8 m/s upstream, midstream and downstream by 2018. The ACD accumulated along the main canal, and the annual average ACDs of downstream were the highest, ranging from 366.17 to 462.95 × 10⁴ cells/L. The joint early-warning method considering both water temperature and flow velocity conditions is an effective way for algal proliferation risk warning management. When water temperatures of the upstream, midstream, and downstream were below 26, 26, and 23°C, respectively, the algal proliferation risk can be controlled under 50% by the flow velocity at 0.3 m/s; otherwise, the flow velocity needs to be regulated higher than 0.8 m/s. In order to keep the midstream and downstream avoid abnormal algal proliferation events (ACD ≥ 500 × 10⁴ cells/L), the corresponding ACDs of the upstream and midstream need to be controlled lower than 319 × 10⁴ cells/L and 470 × 10⁴ cells/L, respectively. This study provides a scientific reference for the long-distance water diversion project’s algal control and environmental protection. The proposed coupling Vine Copula models can also be widely applied to multivariate risk analysis fields.

1. Introduction

With rapid population growth and economic development in recent decades, the water resource conflict between regional water supply and consumption has intensified and become an important global challenge (Cai et al., 2021; Chen et al., 2021). Inter-basin water diversion campaign has become one of the essential engineering measures to alleviate the uneven spatiotemporal distribution of water resources (Zhuang, 2016). Over 160 trans-basin water diversion projects have been built worldwide, for instance, the Middle Route of the South-to-North Water Diversion Project of China (MRSNWDPC), i.e., the world’s longest trans-basin water diversion project (1,432 km), which has delivered a total of more than 60 billion m³ of freshwater resources to northern China since December 2014, directly benefiting more than 100 million people.

Although inter-basin water diversion activities can bring significant economic, ecological, and social benefits (Zhao et al., 2017), such mega hydro projects are often controversial due to the complex anthropogenic imposed on the hydrology, water quality, and ecology of relevant areas (Yang et al., 2018; Dai et al., 2020). The above debate focuses on water quality safety management and the risk of water deterioration before and after project operations. Compared with natural riverine water bodies, despite using the long-distance open canal for water supply, the strict policies of waste discharge prohibition can still eliminate the adverse impact of point source pollution on water bodies, resulting in long-term low nutrient concentrations and low pollution loadings for water quality in the MRSNWDPC (Nong et al., 2020). However, non-point source pollution (e.g., atmospheric dry and wet deposition, endogenous release) may be a major input for nutrient accumulation (e.g., nitrogen and phosphorus) and water quality deterioration in water bodies where long-distance open canals as the primary engineering measure (Zeng et al., 2015; Huang et al., 2019), and the nutrient increase often leads to abnormal algal proliferation (Wasson et al., 2017). On the other hand, the complex algal life cycles (death and decay processes) accompanied by the ecotoxin release can easily become a significant risk factor affecting water quality safety in trans-basin diversion projects (Zhu J. et al., 2021), which has aroused widespread concern and research all over the world (Zhu Y. et al., 2021; Long et al., 2022). Some researchers found that the algal proliferation of the MRSNWDPC may be the main driving factor leading to the chemical oxygen demand increase (Wang C. et al., 2022), causing unacceptable odor and water quality degradation (Wang Y. et al., 2022). However, most of the relevant research is limited to the analysis of algae growth trends and environmental driving factors, and there is still a lack of specific research on how to control algal proliferation.

Abnormal algal proliferation, which marks the precursor of eutrophication, has been widely studied in low-flow velocity, high-depth water column bodies such as lakes and reservoirs (Zou et al., 2020; Li et al., 2021). A series of achievements have been reported regarding key limiting nutrient factors and ecological control measures in algal proliferation (Cheng et al., 2019; Brookfield et al., 2021). For instance, there has become a worldwide consensus that nutrient accumulation (primarily nitrogen and phosphorus) promotes algae growth and reproduction (Huo et al., 2019; Yan et al., 2021). However, early-warning and control strategies of eutrophication in riverine bodies often differ from those in lakes and reservoirs. Many studies pointed out that algae growth could strongly respond to the hydrological condition changes, and the flow velocity can also affect the spatiotemporal distribution of algal cell density (Liang et al., 2013). For instance, regulating flow velocity can lead to sediment resuspension, nutrient redistribution, and spectral composition changes in the water column, thus affecting light energy utilization and effectively inhibiting algal growth (Shen et al., 2021). Low water exchange rates can provide a stable environment for algal growth, thus affecting the migration, diffusion, and accumulation processes (Xin et al., 2019). Additionally, the high-velocity condition can somewhat inhibit the algae growth and even cover up the promotion effects of nutrient loading improvement (Tian et al., 2022). Therefore, it is more important for riverine bodies to comprehensively consider hydraulic regulation and nutrition control strategies to reverse the eutrophication problem efficiently. Since the flow velocity of water bodies in inter-basin water diversion projects can be regulated by sluices (Li et al., 2022), it is feasible to control algal proliferation by changing hydraulic conditions.

Worryingly, the abnormal algal proliferation in the inter-basin water diversion project is more complicated than in rivers. The hydraulic sluice regulation of those water supply projects makes the flow discharge and water level change more frequently than in rivers, speeding up the algal life cycle processes (Zhang et al., 2019). Despite some research claiming that the nutrients in the project’s water bodies were kept at a low level than natural rivers and would not cause algal blooms (Chen et al., 2022), abnormal algal proliferation continues to be observed in the projects (Cattaneo et al., 2013; Li et al., 2014). Previous studies have been limited to the effects of individual key or specific indicators on algae growth, which may hide the factual information about the multivariate risk of abnormal algal proliferation (Xin et al., 2020). The quantification and identification of multi-relationships are necessary to fully explore the real-world ecology-environment-hydrological nexuses (Zhang Y. et al., 2015; Zhao Q. et al., 2019). More deepening studying is needed considering the complex hydrodynamic conditions and low nutrition characteristics of the water environment in those projects (Jing et al., 2022).

Based on the above discussion and gap, the main objectives of this research are: (1) to take the Middle Route of the South-to-North Water Diversion Project of China (MRSNWDPC) as a case and explore the critical factors of algal proliferation in the inter-basin water diversion project under different hydraulic and water quality conditions, (2) to propose a multivariate risk early-warning model and management strategies of algal proliferation under different hydrological and water quality scenarios, (3) to quantify the non-linear effects of algal growth in different locations on the spatial distribution of downstream ecological risks. This study can provide a scientific basis and technical approach for understanding the relationship between planktonic algae and key hydrological and water quality indicators in long-distance water diversion projects, contributing to the adaptive early-warning policies and water quality management strategies for algal bloom risk in different water environments worldwide.

2. Study area and data collection

The Middle Route of the South-to-North Water Diversion Project of China (MRSNWDPC, 32°40′–39°58′ N, 111°42′–116°16′ E) originated from the Danjiangkou Reservoir, the largest artificial freshwater lake in Asia, and northward 1,432 km to Tuancheng Lake, Beijing, the capital of China, since December 2014. As the water quality security of MRSNWDPC is related to the human health of nearly 100 million people, the water quality is monitored by the Construction and Administration Bureau (CAB) of the Middle Route of the South-to-North Water Diversion Project of China as part of the national monitoring program. The CAB has a comprehensive quality assurance/quality control (QA/QC) system to ensure that the data obtained from the national monitoring program is reliable and of sufficient accuracy and precision.

Six representative water quality monitoring stations in the main canal were selected for this study (Figure 1 and Table 1). These stations were selected based on the availability of a full set of water quality indicator data and the coverage of algal cell changes throughout the important cities in the MRSNWDPC. Water samples were collected monthly from January 2016 to December 2018 by the CAB’s professional team according to the environmental quality standards for surface water (China, 2002a). All the sampling date are determined by the CAB and in a unified arrangement. According to the actual weather conditions of each station, the water samples were collected on sunny or cloudy days to minimize the interference of rainfall, strong wind, and runoff on the data and ensured the reliability of the data. All water quality data used in this study were obtained from the CAB with official permission. A data matrix, including algal cell density (ACD, cell/L), pH, dissolved oxygen (DO, mg/L), five-day biochemical oxygen demand (BOD₅, mg/L), ammonia nitrogen (NH₃-N, mg/L), total phosphorus (TP, mg/L), total nitrogen (TN, mg/L), fluoride ( $F^{-}$ , mg/L), water temperature (WT, °C), sulphate ( ${SO}_{4}^{2 -}$ , mg/L), faecal coliform bacteria (F. coli, colonies/L), mercury (Hg, μg/L), arsenic (As, μg/L), copper (Cu, μg/L), zinc (Zn, μg/L), selenium (Se, μg/L), flow velocity (V, m/s), and flow discharge (Q, m³/s), was used for the first step of analysis. Since the water quality of the main canal maintained an excellent status (Nong et al., 2020; Tang et al., 2022), and most of the heavy metal elements and physicochemical indicators had narrow vary ranges, the WT, Q, V, DO, NH₃-N, TP, and TN were selected as the basic indicators that may affect ACD in the main canal based on the pre-review (Song et al., 2018; Zhang et al., 2018; Gao et al., 2022). The ACD was obtained by collecting 1,000 mL of water samples at 0.5 m below the river center’s surface and adding 10–15 mL of Lugol’s iodine solution. The WT and DO were obtained by multi-parameter probe measurements in-situ. For NH₃-N, TP, and TN, water samples were collected at the sampling station and taken in a cool box to the laboratory, where the NH₃-N was determined by the spectrophotometric method with nano reagent, the TP was determined by the potassium persulfate digestion method, and the TN was determined by the potassium persulfate oxidation method (China, 2002b). The monitoring standards and sampling methods can be found in the Ministry of Ecology and Environment of China¹ (in Chinese) and Wu et al. (2018).

FIGURE 1

Figure 1. Location of the monitoring station of the MRSNWDPC in this study [Shahenan (SS); Zhengwan (ZW); Zhanghebei (ZN); Da’anshe (DS); Xiheshan (XS); Huinanzhuang (HN)].

TABLE 1

Table 1. The location of monitoring stations of the MRSNWDPC in this study.

3. Materials and methods

3.1. Vine Copula function

The Copula function is a classical tool to connect multivariate distribution with its univariate probability function (marginal distribution). One of the most important advantages is that it can effectively address the construction issues from marginal to joint distributions for multiple independent variables (Liu et al., 2021). Sklar (1959) first gave this famous theorem as follows.

For a d-dimensional random variable $X_{1}, X_{2}, \dots, X_{d}$ with continuous marginal distribution, thus the joint distribution can be expressed by:

\begin{matrix} F (x_{1} x_{2} \dots x_{d}) = C (F_{1} (x_{1}), F_{2} (x_{2}), \dots, F_{d} (x_{d})) \\ = C (u_{1} u_{2} \dots u_{d}) \end{matrix} (1)

Where $F$ is the corresponding joint cumulative distribution; $F_{i} (x_{i})$ is the marginal cumulative distribution of the random variable $x_{i}$ ; $u_{i} = F_{i} (x_{i})$ , where $i = 1, 2 \dots, d$ . The corresponding density function is (Wan and Li, 2019):

\begin{array}{l} f (x_{1} x_{2} \dots x_{d}) = \prod_{i = 1}^{d} f_{i} (x_{i}) \cdot c (u_{1} u_{2} \dots u_{d}) & (2) \end{array}

Where $f$ is the joint density distribution; $f_{i} (x_{i})$ is the marginal density distribution of the random variable $x_{i}$ , where $i = 1, 2 \dots, d$ ; and $c$ is the density of the Copula function. More details about the Copula can be found in Khozeymehnezhad and Nazeri-Tahroudi (2020) and Nazeri-Tahroudi et al. (2022). Our article mainly introduced and used the six typical Copula functions that were widely applied (Supplementary Table S1). Additionally, there are specific “rotations” structure in some Copula functions, and the corresponding density function transformation can be found in Table 2.

TABLE 2

Table 2. Density function of Copula related to the original morphology after rotation.

Traditional Copula functions are widely used to solve problems in two dimensionals or, at most, three dimensions, and the simulation results need improvement when solving problems in higher dimensions. To address the above issues, Bedford and Cooke (2001) developed the Vine Copula based on the concept of Regular vine (R-vine) as a graphical model of conditional dependence, where Canonical vine (C-vine) and Drawable vine (D-vine) are two unique forms in R-vine Copula (Liebscher, 2008). The principle of Pair-Copula construction (PCC) is to decompose n-dimensional multivariate density functions into $\frac{n (n - 1)}{2}$ binary Copula density functions, thus making it easier to model high-dimensional problems (Kraus and Czado, 2017). Among them, the Vine Copula structure includes three parts: nested tree $T_{i}$ , node $N_{i}$ , and edge $E_{i}$ . For the d-dimensional R-vine Copula (Dong et al., 2022), there are a set of trees $T_{i}$ ( $i = 1, 2, \dots, d - 1$ ) with the structure as follows:

1. T₁ has nodes N₁ = 1, 2, ⋯, d and edges E₁.

2. For i =2, ⋯, d–1, the tree T_i has nodes N_i = E_i–1.

3. Two edges in tree $T_{i}$ are associated to tree $T_{i + 1}$ if they co-use a node in the tree.

A structural diagram of the Vine Copula can be found in Brechmann and Schepsmeier (2013), while the equations of the probability distributions for the Vine Copula are shown in Eq. (3)–(5).

The expression of R-vine Copula can be defined as follows:

\begin{array}{l} f (x_{1} x_{2} \dots x_{n}) = \prod_{k = 1}^{n} f_{k} \prod_{i = 1}^{n - 1} \prod_{e \in E_{i}} c_{j (e), k (e) | D (e)} \\ (F (x_{j (e)} | x_{D (e)}), F (x_{k (e)} | x_{D (e)})) \end{array} (3)

The expression of C-vine Copula can be defined by:

\begin{array}{l} f (x_{1} x_{2} \dots x_{n}) = \prod_{k = 1}^{n} f_{k} (x_{k}) \prod_{j = 1}^{n - 1} \prod_{i = 1}^{n - j} c_{j, j + i | 1, \dots, j - 1} \\ (F (x_{j} | x_{1}, \dots, x_{j - 1}), F (x_{j + i} | x_{1}, \dots, x_{j - 1})) \end{array} (4)

The expression of D-vine Copula can be expressed as:

\begin{array}{l} f (x_{1} x_{2} \dots x_{n}) = \prod_{k = 1}^{n} f_{k} (x_{k}) \times \prod_{j = 1}^{n - 1} \prod_{i = 1}^{n - j} c_{i, i + j | i + 1, \dots, i + j - 1} \\ (F (x_{i} | x_{i + 1}, \dots, x_{i + j - 1}), F (x_{i + j} | x_{i + 1}, \dots, x_{i + j - 1})) \end{array} (5)

Where $f_{k}$ is the edge density function, $\overset{n - 1}{\underset{j = 1}{Π}}$ denotes trees, and $\overset{n - j}{\underset{i = 1}{Π}}$ denotes edges.

The corresponding construction process of the Vine Copula model took the following steps:

1. The Kendell coefficient τ between two variables is calculated, and the structure of Tree 1 in Vine Copula is determined using the maximum spanning tree with τ as the weight.

2. The optimal Pair-Copula among the variables in Tree 1 is determined by the Akaike information criterion (AIC), and the parameter values of the corresponding Copula are calculated using maximum likelihood method.

3. The observed values of Tree 2 are obtained from the parameter values of Tree 1 and Eq. (6), and the optimal Copula of Tree 2 variables is obtained by repeating the previous step.

\begin{array}{l} h (x, v, Θ) = \frac{\partial C_{x v} (x, v, Θ)}{\partial v} & (6) \end{array}

where $Θ$ is the parameter set of the Copula function and $C_{x v}$ denotes the binary Copula distribution function.

4. Repeat steps 2 and 3 until the final layer of the tree is calculated.

5. The Vine Copula parameters were corrected using the maximum likelihood method.

3.2. Latin hypercube sampling

The Latin Hypercube Sampling (LHS) is a method of approximate random sampling from multivariate parametric distributions, which is a stratified sampling technique widely used to generate random samples of parametric values from multi-dimensional distributions (Zhao Y. et al., 2019). The most important advantages of the LHS are the effect of filtering the variance associated with the additive component of the transformation, and the ability to provide uniform solution space coverage and stable output with a smaller number of simulations (Huang et al., 2018). Since the LHS can make full use of the abundant computational resources to simulate a large number of representative samples with strong correlation and cover the whole search space (Feng et al., 2020), a coupling method was proposed to combine the LHS and the Vine Copula model. The model can better handle the dependencies between multiple factors relative to the common linear correlation coefficients (Huang et al., 2020), and compared with machine learning methods such as decision trees, this model will not have the phenomenon of over-fitting and ignoring attribute correlation (Schepsmeier, 2019). In this study, the dependency structure of indicators in each canal section is constructed by the Vine Copula model. Then the data set is expanded by coupling LHS with the Vine Copula model and generating many pseudo-random numbers that meet the high-dimensional dependence of multiple factors, which are correlated with variable correlations in the risk analysis reliably fitted. Further details on LHS are available at (Shields and Zhang, 2016). In this paper, the Vine Copula model $F (X_{1} X_{2} \dots X_{n})$ with n-dimensional variables $(X_{1} X_{2} \dots X_{n})$ was established and an $M \times n$ matrix $(ω_{1} ω_{2} \dots ω_{n})$ follow the uniform distribution of [0,1] generated by the LHS were obtained, where $ω_{i} (i = 1, 2, \dots, n)$ was defined:

{\begin{matrix} ω_{1} = F (X_{1}) \\ \begin{matrix} ω_{2} = F (X_{2} | X_{1}) \\ \begin{matrix} ω_{3} = F (X_{3} | X_{1}, X_{2}) \\ \dots \dots \end{matrix} \end{matrix} \\ ω_{n} = F (X_{n} | X_{1}, X_{2}, \dots, X_{n - 1}) \end{matrix} (7)

Where in each column of data, $U_{i} (i = 1, 2, \dots, n)$ represents the marginal distribution function value of $X_{i} (i = 1, 2, \dots, n)$ . After the LHS simulation, the algae risk of samples was analyzed according to the Bayesian theory (Fan et al., 2020), and the detailed modelling processes are shown in Figure 2.

FIGURE 2

Figure 2. Technical flow chart of algal proliferation risk analysis.

4. Results

4.1. Water quality variation and key water quality factors

4.1.1. Water quality variation characteristic

The statistical summary of the ACD, hydrological factors (Q and V), and water quality indicators (DO, TP, TN, NH₃-N, and WT) of different sections of the MRSNWDPC during the monitoring period can be found in Table 3.

TABLE 3

Table 3. Analysis of water quality parameters in each canal section of the MRSNWDPC from 2016 to 2018 (Avg., average; S.D., standard deviation).

Temporally, the annual average values of the ACD in all sections showed a decreasing trend yearly; spatially, the ACD gradually increased from upstream to downstream along the main canal, and the annual average ACDs of downstream were the highest, ranging from 366.17 to 462.95 × 10⁴ cells/L. There are no significant differences in flow velocity among all sections, which means that although the flow discharge decreased along the canal, the adjustability of the flow velocity can be kept consistent. In addition, the velocity of each section tends to increase year by year, and the maximum velocity will exceed 0.8 m/s by 2018. The WT gradually decreased from upstream to downstream, with a narrow range of 16.6–17.4°C. The annual average concentrations of DO showed a gradual increase from upstream to downstream, while the concentrations of TN have opposite trends. The TP concentration in the upstream is higher than that in the midstream and downstream, and there is no significant spatial difference between the TP concentration in the midstream and downstream. In addition, the average annual concentration of NH₃-N in each canal section is stable, and there is no obvious spatial difference. The lowest TN/TP ratios were in the upstream section (ranged from 25.69 to 44.14), and the middle and downstream (ranged from 58.86 to 101.16, and 64.83 to 97.75, respectively) were significantly higher than those. The nutrient loadings were kept low, with ammonia nitrogen (NH₃-N) ranging from 0.034 to 0.051 mg/L and TP ranging from 0.011 to 0.039 mg/L, respectively. In general, all water quality indicators met the Class I standard, except for TP at the upstream section, which met the Chinese Class II standard, indicating that the MRSNWDPC has excellent water quality conditions and low nutrient loading. In addition, TN (0.778–1.235 mg/L) will not be included in the MRSNWDPC for study according to the Chinese Environmental Quality Standards for Surface Water (China, 2002a).

4.1.2. Relationship between key factors and algal proliferation

The result of the nonlinear relationships between the ACD and the key hydrological and water quality indicators in each canal section based on the multiple regression model can be found in Figure 3, where the R² (coefficient of determination) is calculated by the ordinary least squares (OLS).

FIGURE 3

Figure 3. Response of algal cell density to water quality and hydrological factors variations [Upstream Q-ACD (A), Upstream V-ACD (B), Upstream WT-ACD (C), Upstream DO-ACD (D), Upstream TP-ACD (E), Upstream TN-ACD (F), Upstream NH3-N-ACD (G), Midstream Q-ACD (H), Midstream V-ACD (I), Midstream WT-ACD (J), Midstream DO-ACD (K), Midstream TP-ACD (L), Midstream TN-ACD (M), Midstream NH3-N-ACD (N), Downstream Q-ACD (O), Downstream V-ACD (P), Downstream WT-ACD (Q), Downstream DO-ACD (R), Downstream TP-ACD (S), Downstream TN-ACD (T), Downstream NH3-N-ACD (U)].

The ACD showed a significantly increased trend first and then decreased concerning the increase of Q and V, except for the flow velocity upstream (Q: R² = 0.130 to 0.360, p < 0.05; V: R² = 0.247 to 0.272, p < 0.001). When the WT<13.6°C, the ACD upstream showed a decreased response to the increase of WT, while the ACD of the midstream and downstream showed a simultaneous increase to the WT (R² ≥ 0.480, p < 0.001). The ACD has a negative response to DO concentrations raised both in upstream and downstream (R² = 0.200 and 0.400, p < 0.001), while the ACD first increased and then decreased in the midstream (R² = 0.360, p < 0.001). Except for the TN downstream, the ACD showed no significant response to the variations of TN, TP, and NH₃-N (R² ≤ 0.059, p > 0.05). Based on the above results, the Q, V, WT, and DO were finally determined as the key factors constructing the multivariate risk model with ACD.

4.2. Multivariate risk model construction

The three-dimensional (3D), four-dimensional (4D), and five-dimensional (5D) Vine Copula models of ACD and key water quality and hydrological variables (Q, V, WT, and DO) were constructed, respectively. The data set of the Vine Copula model was simulated by the LHS, and the risk of ACD under different environmental conditions was analysed based on the Bayesian theory. Previous studies have reported that the limited influence of DO on algae proliferation in the MRSNWDPC, and the flow velocity has strong relationship with flow discharge (Zhang et al., 2022). Considering the actual project operation, this study will mainly explore the influence of WT and V on the ACD variations.

4.2.1. Determination of marginal distribution

Nine widely used distribution functions (Normal distribution, Lognormal distribution, Extreme Value distribution, Generalized Extreme Value distribution, Logistic distribution, Pearson III distribution, Weibull distribution, Gamma distribution, and Exponential distribution) were selected to fit the data, and the maximum likelihood method was used to evaluate the parameters. The Kolmogorov–Smirnov (K–S) test was used to test the good fitting of marginal function, and the appropriate marginal distribution was selected by the AIC. Table 4 shows the results of selected marginal functions and the corresponding parameters. It can be seen that both the Generalized Extreme Value distribution and Pearson III distribution were selected 4 times, and the optimal marginal functions of different indicators differed at different locations, indicating the importance of selecting the marginal distribution function by comparison step.

TABLE 4

Table 4. Optimal marginal distribution selection of indicators in each canal section of the MRSNWDPC in this study.

4.2.2. Vine Copula model construction

The optimal Vine Copula models for the upstream, midstream, and downstream of MRSNWDPC in different dimensions are shown in Table 5. The C-vine, D-vine, and R-vine were all applied to construct joint distributions based on the AIC, Bayesian information criterions (BIC), and log-likelihood (L-lik). It can be seen in Table 5 that the optimal Vine Copula models of all canal sections selected the R-vine structure. The Vine Copula structure with the lowest AIC and BIC and highest L-lik can be selected as the optimal model (Tahroudi et al., 2022), details can be found in Supplementary Tables S2–S4, and the corresponding Vine Copula structure can be found in Supplementary Figures S1–S3.

TABLE 5

Table 5. Optimal Vine Copula model selection in each canal section of the MRSNWDPC in this study.

4.2.3. Risk analysis of abnormal algal proliferation

Based on the optimal Vine Copula structure and LHS obtained in different dimensions, a set of data sets $[R_{A C D}, R_{W T}, R_{V}]$ , $[R_{A C D}, R_{W T}, R_{V}, R_{DO}]$ , and $[R_{A C D}, R_{W T}, R_{V}, R_{DO}, R_{Q}]$ with the size of 1,000,000 × M (M = 3, 4, and 5) were simulated for each canal section, respectively. The risk probability corresponding to ACD variation for each section under different WT and V conditions was calculated according to Eq. (8).

\begin{array}{l} F (A C D > x | t_{i} < W T < t_{i + 1}, v_{j} < V < v_{j + 1}) \\ = \frac{F (A C D > x, t_{i} < W T < t_{i + 1}, v_{j} < V < v_{j + 1})}{F (t_{i} < W T < t_{i + 1}, v_{j} < V < v_{j + 1})} = \frac{n}{N} \end{array} (8)

Where $i = 1, 2 \dots, 6; j = 1, 2 \dots, 5; N$ denotes the number of samples of $t_{i} < W T < t_{i + 1}$ and $v_{j} < V < v_{j + 1}$ ; n denotes the number of samples satisfying both $t_{i} < W T < t_{i + 1}, v_{j} < V < v_{j + 1}$ and $A C D > x .$

It can be seen from Supplementary Figures S4–S9 that when the water temperature is within a specific interval, the probability of ACD proliferation of the upstream and midstream will decrease with the increased flow velocity. For instance, when the water temperature of upstream between 26 and 29°C, the ACD had a 61.4% probability of exceeding 300 × 10⁴ cells/L under the V under 0.3 to 0.4 m/s condition, when the V raised to 0.7 to 0.8 m/s, the corresponding rate reduced to 30.9%. However, in the downstream section, the ACD variations showed weak relationships with flow velocity changes. It can be observed that the flow velocity changed from 0.3 to 0.6 m/s had almost no effect on the probability of exceeding a specific ACD value under the same water temperature condition.

4.3. Spatial relationships of ACDs in different canal sections

The histogram, scatter plot, and Kendall coefficient τ after normalizing the ACD data of each channel section were presented in Figure 4. It can be seen that the ACDs did not show normal distributions in each section. There were strong correlations between the ACD of adjacent sections, with the Kendall coefficient τ (p ≤ 0.05) of “up to midstream” and “mid to downstream” were 0.68 and 0.87, respectively.

FIGURE 4

Figure 4. Spatial relationships of ACDs in different canal sections (*p ≤ 0.05; **p ≤ 0.01; ***p ≤ 0.001).

The responses of the ACD between different canal sections were analyzed based on binary Copula models (Supplementary Table S5), and the conditional probabilities of downstream ACD variation under the upstream ACD change were calculated based on the Bayesian theory according to Eq. (9).

\begin{array}{l} F (AC D_{i} > x | y_{1} \leq AC D_{j} < y_{2}) = 1 - \frac{C (x, y_{1}) - C (x, y_{2})}{F_{AC D_{j}} (y_{2}) - F_{AC D_{j}} (y_{1})} & (9) \end{array}

Where i = Mid, Down; corresponding to j = Up, Mid.

Considering the ACD variation range of each section, the conditional probability of ACD_Mid (ACD_Down) changes under the different ACD_Up (ACD_Mid) is shown in Figure 5. With the ACD increased in upstream, the risk probability of ACD_Mid increased simultaneously, where the ACD_Mid only had a 2.48% probability exceeding of 500 × 10⁴ cells/L when the ACD_Up ≤ 50 × 10⁴ cells/L. However, when the ${A C D}_{Up}$ raised to 350–400 × 10⁴ cells/L, that risk probability was higher than 62%. Similar phenomena can also be observed in the mid to downstream section. From the spatial prospect, the ACD control using the flow velocity regulation strategy should be carried out in the upstream and middle sections, which can achieve efficiency.

FIGURE 5

Figure 5. The conditional probability of ACD variations on adjacent canal sections [(A,B) are defined as the effects of ACD_Up on ACD_Mid and ACD_Mid on ACD_Down, respectively].

5. Discussion

5.1. Key driving factors of ACD variations

During the study period, the ACD in MRSNWDPC was maintained at the level of 10⁶ cells/L, in which the dominant algae composition was of Chlorophyta, Cyanophyta, and Bacillariophyta (diatoms). Chlorophyta and Cyanophyta were the dominant species in summer, and the Cyanophyta and Bacillariophyta were the dominant species in the rest of the seasons (Zhang et al., 2021c). Previous studies have shown that water temperature is one of the critical factors that can stimulate algae growth (Zhao et al., 2018; Mei et al., 2022), and the growth rate of Chlorophyta and Cyanophyta can reach the peak when surface water temperature within 25–35°C (Huisman et al., 2018). This study also showed similar results (Figure 6), reporting that the high risk of algal proliferation happened when the upstream, midstream, and downstream had the WT higher than 26, 26, and 23°C, respectively. Those phenomena can be attributed to the dominant species Chlorophyta and Cyanophyta which can have better reproductive ability in higher water temperature environments (Li R. et al., 2015). Considering the seasonal characteristic of the MRSNWDPC (Zhang et al., 2021a), the prevention and control management of ACD in the canal should be strengthened in summer and autumn.

FIGURE 6

Figure 6. Early-warning thresholds for algal proliferation control under different conditions at upstream (A), midstream (B), and downstream (C).

Flow velocity is an important factor affecting algal habitats (Li Z. et al., 2015). Under suitable velocity conditions, the algae growth rate will be promoted, while high flow velocity can also inhibit algal growth (Zhou et al., 2018). The flow conditions change causes the sediment to resuspend, leading to the nutrient re-distribution and spectrum composition variation. The above processes affect light utilization and the nutrient absorption of algae growth. Additionally, under high-velocity conditions, some algae cells will be damaged by the impact of water flow and accelerate to extinction process (Zhang H. et al., 2015). For the upstream section of the MRSNWDPC, the ACD proliferation risk can be reduced by regulating the flow velocity above 0.6 m/s when the water temperature was higher than 26°C, which may be due to the high flow velocity inhibiting the growth of algae cells and nutrient uptake, and reducing the photosynthetic pigment content and PSII (photosystem II) activity (Li et al., 2022). However, it can be seen that the flow velocity changes between 0.3–0.7 m/s in the midstream and downstream have little effect on the ACD variations (Supplementary Figures S7–S12), which may be due to the increase of algae in the midstream and downstream, which makes a large number of algae adhere to the bend of the canal and the side walls (Hu et al., 2022), resulting in the regulation of algae by relatively low flow velocity showed no significant effects. Therefore, the flow velocity needs to be pre-regulated for the ACD control before the WT exceeding 26 and 23°C, respectively. When the water temperature exceeds the warning threshold, increasing the flow velocity cannot achieve satisfactory results in the middle and downstream sections.

5.2. Algal-controlled threshold for different sections

The early-warning thresholds of algal cell density in the upstream, midstream, and downstream sections of the MRSNWDPC were set to 300 × 10⁴ cells/L, 500 × 10⁴ cells/L, and 500 × 10⁴ cells/L, respectively, according to the actual situation of the project and previous study (Zhang et al., 2022). The event in which the ACD has a 50% probability of exceeding the threshold was defined as “algal proliferation risk.” In this research, it can be seen that the ACD in the relative upstream location can significantly affect the ACD of the downside. When the ACD_Up < 182 × 10⁴ cells/L, the ACD_Mid has a 20% probability exceeding 500 × 10⁴ cells/L. According to the binary Copula model calculation, the ACD_Up = 319 × 10⁴ cells/L and ACD_Mid = 470 × 10⁴ cells/L were selected as the early-warning thresholds for the upstream and midstream section, respectively, which were about ±6% differences from previous study (300 × 10⁴ cells/L and 500 × 10⁴ cells/L, respectively), reflecting that the algal controlled threshold for ACD_Up and ACD_Mid determined by the Copula model was appropriate. For the practical control and management of the project, the early-warning thresholds of ACD in the upstream, midstream, and downstream sections of the MRSNWDPC can be proved reasonably set to 300 × 10⁴ cells/L, 500 × 10⁴ cells/L, and 500 × 10⁴ cells/L, respectively. Additionally, it was found that the average concentration of ACD in summer often exceeded 500 × 10⁴ cells/L in downstream (Nong et al., 2021). Therefore, relevant management can be conducted by adjusting flow velocity higher than 0.8 m/s in the upstream section before summer (Atazadeh et al., 2021).

5.3. Model uncertainty

This research introduced the Vine Copula model to construct multivariate joint distributions for algal proliferation risk analysis. The effects of water temperature and flow velocity on the ACD variations in different dimensions were analyzed. There were significant differences between the results of the binary Copula models (Zhang et al., 2021b) and this study, which may be due to the neglect of the multivariate correlations in the previous research. Therefore, multidimensional joint risk modelling should be conducted by considering the comprehensive effect of multiple factors simultaneously to improve the realism of the simulation results. It was found that the fitting results of R-vine Copula were better than that of C-vine and D-vine Copula, which may be attributed to the fact that the R-vine Copula function can better maintain the statistical characteristics of Kendall and Spearman correlation coefficients of the original sequence (Zeng et al., 2022). The proposed coupling method using the LHS and the Vine Copula model showed robustness and comparable ability. The simulation results were consistent with the Monte Carlo simulation by Zhang et al. (2022). However, when analyzing the probability under special conditions such as the ACD > 0 cells/L (probability = 100%), the results of this study were more in line with the actual situation. Additionally, the risk intervals of the ACDs showed little difference compared with the Monte Carlo simulation (Zhang et al., 2022), which can be attributed to the uniform sampling characteristic of the LHS method (Zhang W. et al., 2015). To sum up, the proposed coupling method in this research can effectively reduce sampling and simulation uncertainty and can be generalized and applied to more multi-dependent structure analysis and risk assessment.

6. Conclusion

This study conducted a spatiotemporal algal proliferation risk analysis in the MRSNWDPC, considering the effects of hydrological and water quality indicators. A comprehensive coupling method based on the Vine Copula model and Latin hypercube sampling (LHS) technique was proposed. The algal proliferation risk under different environmental conditions in the main canal was analyzed. The main conclusions are:

1. All the water quality indicators during the monitoring periods from 2016 to 2018 meet the Class I or II standard, showing the excellent and low nutrient loadings of water quality status of the South-to-North Water Diversion Project of China. However, the algal cell density gradually increased from upstream to downstream along the main canal, and the annual average ACDs in downstream were the highest.

2. The joint early-warning system considering water temperature and flow velocity conditions is an effective way of algal proliferation control. The prevention and control management of ACD in the canal should be strengthened in summer and autumn. From the spatial prospect, the ACD control using the flow velocity regulation strategy should be carried out in the upstream and middle sections, which can achieve efficiency.

3. When water temperatures of the upstream, midstream, and downstream were below 26, 26, and 23°C, respectively, the algal proliferation risk can be controlled under 50% by the flow velocity at 0.3 m/s; otherwise, the flow velocity needs to be regulated higher than 0.8 m/s, which can improve the water quality of the main canal while controlling algae via flow velocity regulation. In order to keep the midstream and downstream avoid abnormal algal proliferation events (ACD ≥ 500 × 10⁴ cells/L), the corresponding ACDs of the upstream and midstream need to be controlled lower than 319 × 10⁴ cells/L and 470 × 10⁴ cells/L, respectively.

4. The coupling method using the Latin hypercube sampling model and the Vine Copula model can effectively reduce sampling and simulation uncertainty in multivariate risk analysis, and can be generalized and applied to more multi-dependent structure analysis and risk assessment.

This study can provide a practical abnormal algal proliferation management strategy in the MRSNWDPC, which can effectively control the algal growth under different water temperatures by regulating flow in specific canal sections, and the early warning threshold for algal prevention can provide a comparable reference for similar water bodies. The proposed coupling method with the Vine Copula model and the LHS method can also be applied for the multivariate risk analysis field. However, due to limitation of data availability, this study may ignore the impact of some meteorological factors (e.g., rainfall and solar radiation) on algal proliferation. Considering that those projects will continue to supply water for large-scale cities for many years, long-term meteorological, hydrological, water quality, and ecology monitoring is still strongly needed, as well as in-depth multivariate risk analysis in future work.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary material, further inquiries can be directed to the corresponding author.

Author contributions

XN: data collection, investigation, conceptualization, methodology, formal analysis, writing – original draft, and funding acquisition. JZ and YM: software and writing – original draft. LC: formal analysis, draft review and editing, and funding acquisition. CZ: data collection, investigation, and draft review and editing. KB: draft review and editing, and supervision. LCC: formal analysis, draft review and editing, and supervision. All authors contributed to the article and approved the submitted version.

Funding

This research was funded by the Specific Research Project of Guangxi for Research Bases and Talents (No. AD22035185), the Youth Science Foundation of Guangxi (No. 2023GXNSFBA026296), The Belt and Road Special Foundation of the National Key Laboratory of Water Disaster Prevention (No. 2022nkms06), and the National Natural Science Foundation of China (No. 51669003).

Acknowledgments

The authors would like to thank the editors and the reviewers for their valuable suggestions and contributions which significantly improved this article, and would also need to acknowledge the Construction and Administration Bureau of the Middle Route of the South-to-North Water Diversion Project of China that supported the data collection.

Conflict of interest

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

Publisher’s note

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

Supplementary material

The Supplementary material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fevo.2023.1193163/full#supplementary-material

Footnotes

1. ^ https://www.mee.gov.cn/ywgz/fgbz/bz/

References

Atazadeh, E., Gell, P., Mills, K., Barton, A., and Newall, P. (2021). Community structure and ecological responses to hydrological changes in benthic algal assemblages in a regulated river: application of algal metrics and multivariate techniques in river management. Environ. Sci. Pollut. Res. 28, 39805–39825. doi: 10.1007/s11356-021-13546-w