The Mann-Kendall (MK) test (Mann, 1945; Kendall, 1975) and Pettitt (PT) test (Pettitt, 1979) are widely used as non-parametric methods for trend and change-point detection in hydro-meteorological series, respectively.

In order to accurately assess the reservoir impacts on the downstream hydrological variables, in this study, a general rainfall-reservoir index (GRRI) is developed by coupling effective rainfall that drives hydrological variables and reservoir regulation indicator, and it is defined as follows G R R I = 0 , R R = 0 1 − F E R 1 / R R − 1 , 0 < R R ≤ 1 R R , R R > 1 (7) where, F _ER is the probability distribution of effective rainfall, which can be the probability distribution of univariate rainfall event or the joint probability distribution of multivariate rainfall event; RR is an indicator that symbolizes the impact of reservoir regulation, e.g., reservoir index (RI) and sediment trapping efficiency (TE). The relationship between GRRI and RR and F _ER is shown in Figure 2.

Since F E R ∼ U 0,1 , 1 − F E R ∼ U 0,1 , 1 − F E R = G R R I R R / 1 − R R , when 0 < R R ≤ 1 , the cumulative probability distribution function and probability density function of GRRI are F G R R I y G R R I = 0 , y G R R I = 0 y G R R I R R / 1 − R R , 0 < y G R R I ≤ 1 1 , y G R R I > 1 (8) f G R R I y G R R I = R R 1 − R R ⋅ y G R R I R R / 1 − R R − 1 , 0 < y G R R I ≤ 1 0 , y G R R I = 0   o r   y G R R I > 1 (9)

When 0 < R R ≤ 1 , the expectation of the GRRI is E G R R I = ∫ − ∞ + ∞ f G R R I y G R R I ⋅ y G R R I d y G R R I = R R .

For various hydrological variables, different effective rainfall and reservoir regulation indicators are selected to calculate GRRI as covariates for nonstationary hydrological frequency analysis. In this study, the two hydrological variables of flood and annual sediment load are taken as examples, their corresponding GRRIs are calculated according to Eq. 7, which are marked as GRRI_Q and GRRI_S respectively. Specific details on how to identify the effective rainfall and reservoir regulation indicator for flood and annual sediment load are described below.

The covariate-based nonstationary frequency analysis has been widely concerned by hydrologists (Strupczewski et al., 2001; Rigby and Stasinopoulos, 2005; Villarini et al., 2010; Salas and Obeysekera, 2014; Jiang et al., 2015; Xiong et al., 2020; Li et al., 2022), in which distribution parameters are expressed as the function of covariates/explanatory variables. In this study, Weibull (WEI), Gumbel (GU), Gamma (GA), Logistic (LO), Normal (NO), Lognormal (LNO), Generalized extreme value (GEV), and Pearson type III (P-III) distribution are used as candidate distributions to fit hydrological series, the location parameter ( μ t ) and scale parameter ( σ t ) of the distribution are considered to change with the covariates. For the three-parameter GEV and P-III distribution, the shape parameter ( ξ ) is generally regarded as constant in order to reduce the uncertainty caused by parameter estimation due to the sensitivity of ξ . The nonstationary probability distribution models are constructed based on Generalized Addictive Models in Location, Scale, and Shape (GAMLSS) model by using the reservoir regulation indicator (RR) and the general rainfall-reservoir index (GRRI) as the covariate respectively. And the linear and exponential functions are considered to describe the relationship between distribution parameters and covariates. Taking the case that both location parameter and scale parameter change with covariates as an example, the formula of distribution parameters is expressed as follows L i n e a r : μ t = α 0 + α 1 ⋅ R R σ t = β 0 + β 1 ⋅ R R   o r   μ t = α 0 + α 1 ⋅ G R R I σ t = β 0 + β 1 ⋅ G R R I E x p o n e n t i a l : μ t = exp α 0 + α 1 ⋅ R R σ t = exp β 0 + β 1 ⋅ R R   o r   μ t = exp α 0 + α 1 ⋅ G R R I σ t = exp β 0 + β 1 ⋅ G R R I (13) where α 0 , α 1 , β 0   a n d   β 1 are model parameters estimated by the maximum likelihood estimate method. The goodness-of-fit test of probability distribution models is examined by Kolmogorov–Smirnov (KS) test (Massey, 1951), if the p-value of KS test is greater than the significance level of 0.1, it means that the model fits well. The relative fitting qualities of the probability distribution models are evaluated by the Schwarz Bayesian Criterion (SBC; Schwarz, 1978), and the model featured with the smaller SBC value is considered better.

The Yangtze River is the largest river in China measuring approximately 6,397 km in length. It spans the three steps of China’s terrain from west to east, and flows through 11 provinces. The total basin area is about 1.8 million km², accounting for 18.8% of China’s land area. The vast area and complex topography have created the diverse monsoon climate characteristics of the Yangtze River basin. The average annual precipitation is about 1,100 mm, with large inter-annual variation and uneven distribution within the year. The annual runoff of the Yangtze River basin reaches 960 billion m³ and is rich in hydropower resources. The degree of hydropower development and water resource utilization in the basin is relatively high, forming the world’s largest reservoir group with the Three Gorges Reservoir as the core (Xu and Zhang, 2018; Li H. et al., 2020). The basin has a well-developed river network and numerous tributaries, with 8 basins exceeding 80,000 km², 49 basins more than 10,000 km². In this paper, the Min River basin, Wu River basin, Han River basin and the Cuntan station on the main stream (Figure 3) are taken as typical study areas to carry out the nonstationary hydrological frequency analysis.

The Min River is situated 99°37′E-104°38′E and 28°13′N-33°38′N, covering a total area of 135,387 km². It is located in the upper reaches of the Yangtze River with a length of 735 km. The annual average temperature is 15°C and mean annual precipitation ranges from 800 mm to 1,100 mm (Tang et al., 2013). The Wu River, with the coordinates of 104°18′E-109°22′E and 26°07′N-30°22′N, and a catchment area of 87,920 km², is the largest tributary in the upper reaches of the Yangtze River measuring approximately 1,050 km in length. The annual average temperature varies from 13°C to 18°C and mean annual precipitation is from 900 mm to 1,400 mm (Wu et al., 2018; Li R. et al., 2020). The Han River is situated 106°00′E-114°00′E and 30°30′N-34°30′N, covering a total area of 159,000 km². It is located in the middle reaches of the Yangtze River with a length of 1,532 km and is the largest tributary of the Yangtze River. The annual average temperature is 14°C–16°C and mean annual precipitation varies from 700 mm to 1,100 mm (Xiong et al., 2019; Xie et al., 2021). The Cuntan hydrological station is a crucial control station for the main stream of the upper Yangtze River, controlling nearly one-half of the Yangtze River basin (catchment area: 866,559 km²), the annual average temperature is 13°C and mean annual precipitation is 1,004 mm. Many reservoirs in these basins have been built and put into operation, among which, the nonstationary analysis of Cuntan station only considers the large reservoirs with a total storage capacity greater than 1 billion m³, and the information of all reservoirs is listed in Table 1.

In this study, the collected data include streamflow data, sediment data, rainfall data and reservoir data in the study area. Specifically, daily discharge, daily rainfall, and annual sediment load records during 1960–2019 of Gaochang (GC) station in Min River; daily discharge, daily rainfall, and annual sediment load records of Wulong (WL) station in Wu River from 1952 to 2019; daily discharge and daily rainfall records (1956–2019) of Ankang (AK) station and Huangjiagang (HJG) station in Han River; daily discharge, daily rainfall, and annual sediment load records during 1953–2019 of Huangzhuang (HZ) station in Han River; and daily discharge, daily rainfall, and annual sediment load records of Cuntan (CT) station on the main stream of Yangtze River from 1960 to 2019. The streamflow data and sediment data are provided by the Hydrology Bureau of the Changjiang Water Resources Commission, China (http://www.cjh.com.cn/en/), and the rainfall data is obtained from the National Climate Center of the China Meteorological Administration (http://data.cma.cn/site/index.html).

The Mann-Kendall test and Pettitt test are used for trend and change-point testing of annual maximum daily discharge (Q) and annual sediment load (S) of GC, WL, AK, HJG, HZ, and CT station, the results are shown in Table 2. From the results in the table, it can be seen that, except for Q at CT station, the other hydrological variables of the six stations have significant decreasing trends or change-points at the significance level of 0.1. It is not difficult to find that the change-points in the hydrological series are closely linked with the time of the construction of the large reservoirs in the basin. Since the nonstationary diagnosis results of flood series at CT station are not significant, the nonstationary flood frequency analysis of CT station will not be carried out below.

The inter-annual time series of Q and S at GC, WL, AK, HJG, HZ, CT station are displayed in Figure 4. It can be seen from the figure that all hydrological series at the six stations presented downward trends after 1980 or 1990 when large reservoirs were built in the basin, which is consistent with the results of Mann-Kendall test. As displayed in Table 3, both the mean and coefficient of variance (Cv) of each hydrological series change greatly after the change-point compared with that before the change-point. In addition to the mean of Q at CT station reduces by 10.6%, the mean of Q at other stations decreases by 22.8%–60.6%, and the mean of S drops sharply by 47.7%–89.5%.

The stationary hydrological frequency analysis and nonstationary hydrological frequency analysis with RR (RI and TE) or GRRI (GRRI_Q and GRRI_S) as the covariate are carried out separately for Q and S at each hydrological station to investigate how the reservoirs affect downstream flood and sediment load. The summary of the results for fitting the stationary and nonstationary probability distribution models to the Q and S is shown in Table 5. The relative size of SBC values indicates that the nonstationary models are better than the stationary models for Q and S at all stations, and the nonstationary probability distribution models with GRRI as the covariate have better fitting effects than nonstationary models with RR as the covariate, and the results of KS test indicating that the selected models are all reasonable. Therefore, the nonstationary model with GRRI as the covariate is the best model, which is obviously more suitable for the hydrological frequency analysis at GC, WL, AK, HJG, HZ, and CT station. For the nonstationary probability distributions, the location parameter μ t decreases with the increase of RR or GRRI, revealing the decreasing degree of the frequency and magnitude of downstream flood and sediment load due to the reservoir effects.

The quantile curve plots (Figure 8 for flood and Figure 9 for annual sediment load) demonstrate that the magnitude of Q and S decreases since 1980 or 1990 when a large number of reservoirs are constructed upstream of the hydrological station, and hydrological series is well fitted by the best models. Undoubtedly, with the incorporation of the impacts of effective rainfall, the GRRI can more conveniently and accurately capture the occurrence of nonstationarity in the downstream hydrological frequency. Taking the HJG station and WL station as an example, when the Danjiangkou Reservoir, 6.19 km upstream of the HJG, was built in 1967, flood discharge of HJG began to decrease. As more reservoirs were put into operation, flood discharge continued to decrease. However, some relatively large events still occurred, such as 20,700 m³/s in 1975. When the Wujiangdu Reservoir was completed in 1983, sediment load of WL began to decrease. As more reservoirs were constructed, the magnitude of sediment load kept decreasing. Whereas, some relatively large events still occurred, such as 3,575 × 10⁴ t in 1996. Obviously, the occurrence of this phenomenon can be well explained by GRRI, while RR is powerless in this regard.