Uncertainty Quantification of Rainfall-runoff Simulations Using the Copula-based Bayesian Processor: Impacts of Seasonality, Copula Selection and Correlation Coefficient

The outputs of Rainfall-runoff models are inherently uncertain and quantifying the associated uncertainty is crucial for water resources management activities. This study presents the uncertainty quantification of rainfall-runoff simulations using the copula-based Bayesian processor (CBP) in Danjiangkou Reservoir basin, China. The seasonality of uncertainty in rainfall-runoff modeling is explored, and impacts of copula selection and correlation coefficient on uncertainty quantification results are investigated. Results show that the overall performance of the CBP is satisfactory, which provides a useful tool for estimating the uncertainty of rainfall-runoff simulations. It is also demonstrated that the dry season has higher reliability and greater resolution compared with wet season, which illustrates that the CBP captures the actual uncertainty of rainfall-runoff simulations more accurately in dry season. Moreover, the performance the CBP highly depends on the selected Copula function and considered Kendall tau correlation coefficient. As a result, great attention should be paid to selecting the appropriate Copula function and effectively capturing the actual dependence between observed and simulated flows in the CBP-based uncertainty quantification of rainfall-runoff simulations practice.


Introduction
Rainfall-runoff models are important tools for hydrological forecasting, climate change impact assessment and other water-related activities (Liu et al. 2018;Yao et al. 2021). However, the outputs of the rainfall-runoff model are inevitably subject to errors and inherently uncertain resulting from various sources including input, model parameters and structure (Ragab et al. 2020;Moges et al. 2021). Therefore, understanding and quantification of modeling uncertainty is crucial for risk-based decision-making and has received increasing attention from researchers and practitioners (Biondi et al. 2010;Pappenberger et al. 2015;Yang et al. 2020;Zhou et al. 2020).
Numerous approaches for quantifying the uncertainty in rainfall-runoff simulations have been proposed in the literatures, and these approaches are generally fall into four categories according to uncertainty category concerned (Moges et al. 2021): (1) input uncertainty Fraga et al. 2019), (2) model parameters uncertainty (Gopalan et al. 2019;Ghaith and Li 2020), (3) model structure uncertainty (Knoben et al. 2020;Dakhlaoui and Djebbi 2021) and (4) integrated uncertainty (Liu et al. 2018;Zhou et al. 2020). Actually, there is usually interaction among the uncertainties of model inputs, model parameters and model structure, which jointly determine the uncertainty of model outputs (Moges et al. 2021;Zhou et al. 2021). Many hydrologists have highlighted the importance of integrated uncertainty assessment in rainfall-runoff simulations for practical use and devoted their efforts to develop various approaches (Krzysztofowicz and Kelly 2000;Liu et al. 2018;Yao et al. 2019;Zhou et al. 2020).
To the authors' knowledge, Krzysztofowicz and Kelly (2000) developed the Hydrologic Uncertainty Processor (HUP) which becomes the most commonly used approach to assess the integrated uncertainty of rainfall-runoff simulations (Yang et al. 2020). The HUP has been used to various basins and achieved good results after verifying (Biondi et al. 2010;Biondi and De Luca 2013;Han et al. 2019). However, due to the non-Gaussian characteristic of many hydrological variables, the HUP have to rely on Normal Quantile Transform (NQT) and linear-normal hypothesis, which restricts its application scope and convenience. To overcome this limitation, Liu et al. (2018) propose a copula-based Bayesian processor (CBP) where the prior density and likelihood function were explicitly expressed in the original non-Gaussian space. Furthermore, the CBP allows for Gaussian or non-Gaussian marginal distribution of hydrological variables, and a nonlinear and heteroscedastic dependence structure, which has been tested and proven to be reliable and skillful, thus will provide a promising tool for estimating the uncertainty of rainfall-runoff simulations.
Unfortunately, there are some issues still remain unsolved in the present uncertainty quantification of rainfall-runoff simulations using the HUP and CBP. On the one hand, only one set of parameters was used to construct the uncertainty assessment model throughout the year, without considering the variation of parameters with seasons. On the other hand, although some studies (Li et al. 2017;Feng et al. 2019) have explored how the marginal distribution type impacts on uncertainty assessment, the impact of different joint distributions, including copula selection and different correlation coefficients has received less attention in previous studies.
The aims and objectives of this study are: (i) presenting uncertainty quantification of rainfall-runoff simulations using the CBP by applying to a case study in Danjiangkou Reservoir basin, China; (ii) exploring the seasonality of uncertainty in rainfall-runoff modeling; (iii) investigating impacts of copula selection and correlation coefficient on uncertainty quantification results. The remainder of the paper is organized as follows. The methodology is outlined and presented in Sect. 2. Section 3 introduces the case study. Then, the results and analysis obtained are included in Sect. 4. The discussions are provided in Sect. 5. Finally, Sect. 6 draws the main conclusions of the work.

Marginal Distributions of Observed and Simulated Flows
The probability distributions of daily observed and simulated flows involved in our study refer to the flow duration curve, and the daily streamflow is assumed to be a random variable (Castellarin et al. 2004). The Gamma and Log-normal are commonly used distributions in hydrology to fit the daily streamflow series in the literatures (Xiong et al. 2015), thus will be selected as two candidate theoretical probability distributions of observed flow H and simulated flow S in this study. The probability density function of the Gamma and Log-normal are as follows.
where α and β represent the parameters of the Gamma distribution.
where Y and Y represent the parameters of the Log-normal distribution.
The L-moment method (Hosking 1990) and Kolmogorov-Smirnov test are used to estimate parameters and verify marginal distributions. For more details please refer to Liu et al. (2016).

Joint Distributions of Observed and Simulated Flows
The dependence structure between observed flow H and simulated flow S is modelled with a copula function. In this study we concentrate on the bivariate case. By virtue of the copula function, the joint distribution function of observed flow H and simulated flow S can be expressed as (Razmi et al. 2022) where u = G(h) and v = F(s) are marginal distributions of observed flow H and simulated flow S, respectively; θ represents the parameter of a bivariate copula.
The Gumbel-Hougaard (G-H) copula, Frank copula and Clayton copula are three widely used bivariate Archimedean copulas in hydrology, whose mathematical expressions are described in Table 1 (Wu et al. 2021;Li et al. 2022), are adopted as candidates. The Kendall correlation coefficient tau (τ) inverse method and Cramér-von Mises test are adopted to estimate the parameter θ and verify of the bivariate copula. For more details please refer to Dung et al. (2015).

Copula-based Bayesian Processor
Given the deterministic simulated flow S = s , the corresponding Bayesian posterior density function of observed flow H is expressed as in Eq. (4): The prior density used in the CBP is directly the probability density function of H which expressed as g(h) . It can be estimated based on historical data of observed flows using the procedure in Sect. 2.1 in the original space. Given S = s , the likelihood function using the copula function is defined as (Motevali et al. 2021) where With the help of a copula, the (h|s) is immediately expressed as The numerical solution of the posterior density function (h|s) can be calculated by using the Monte Carlo sampling approach. The corresponding posterior distribution function Φ(h|s) is defined by

Performance Criteria Used
The Nash-Sutcliffe efficiency (NSE) and Relative Error (RE) are used as performance metric for deterministic simulations (Al-Juboori 2022), while quantile-quantile (QQ) plot, reliability (α-index), resolution (π-index) and continuous rank probability score (CRPS) for probabilistic simulations (Pappenberger et al. 2015;Liu et al. 2018). For the calculation method and evaluation criteria for each metric, please also refer to the literature mentioned above.

Case Study
The Danjiangkou Reservoir located at the upper Hanjiang River (Fig. 1), the longest tributary of the Yangtze River, is selected as case study. The catchment area is 95,200km 2 and the average annual runoff amount is 40.85 billion m 3 . It is a multipurpose reservoir with main functions, including flood control, water supply, hydropower generation, and irrigation . Hence, accurate simulation of its flows is of great importance. In this study, the variable infiltration capacity (VIC) model (Shen et al. 2020;Sepúlveda et al. 2022) was established to simulate daily flows of Danjiangkou Reservoir basin. The used dataset and parameter calibration procedures of the VIC model are similar with Guo et al. (2009) and the readers can refer to this paper for more information.

Calibration of CBP
Marginal distributions for four random variables (observed and simulated flow in two seasons) needed to be estimated. The estimated parameters of Gamma and Log-normal distributions were presented in Table 2. The statistic D values of K-S tests were also listed in Table 2. It is indicated that the Log-normal distribution has passed the K-S test at the 5% significance level, while rejected for Gamma distribution. For the observed and simulated flow in two seasons, Log-normal distributions are chosen as theoretical fitting distributions.
The τ values of observed and simulated flows are 0.760 and 0.658 in wet and dry season, respectively. The results of the copula parameter θ and Cramér-von Mises test were summarized in Table 3. These results demonstrate that the G-H and Frank copulas have passed. Furthermore, G-H copula was finally selected owing to its better performance (lower S n value) for each season.

Performance Assessment of the CBP
In this study, for a given significance level 0.1, the 90% confidence interval of observed flow H with the 5% quantile and 95% quantile can be calculated. As typical examples, the observed and median discharges, and 90% simulation intervals in wet and dry seasons estimated by CBP in 1983 (wet year), 1981 (normal year) and 1988 (dry year) were presented  Fig. 2. It is revealed that the observed flows are basically included in the 90% simulation intervals, which indicates that the derived 90% simulation intervals can convincingly assess the simulation uncertainty and supply risk information for water-related decision-making.  Figure 3 shows the simulative QQ plots in wet and dry seasons. The result is encouraging and the general performances of simulative QQ plots in each season are pretty good. However, the behavior is found to be different in wet and dry seasons. The CBP somewhat underestimates simulative uncertainty in wet season, while overestimates simulative uncertainty in dry season. The QQ plot for dry season is a bit closer to the 1:1 line than wet season, which demonstrates the CBP performs slightly better in dry season in terms of reliability.
The results of α-index, π-index and CRPS were listed in Table 4. It is clearly shown that both the α-index and π-index values in dry season are larger than those of the wet season. Compared with wet season, the dry season has higher reliability and greater resolution, which illustrates that the CBP captures the actual uncertainty of rainfall-runoff simulations more accurately in dry season. In terms of the CRPS, both wet and dry seasons have lower values than those of deterministic simulations and indicate the effectiveness of probabilistic simulations. The CRPS value in wet and dry seasons is decreased by 34.6% and 43.3%, respectively. All these results reveal that in this case study the performance of the CBP in dry season is superior to the wet season and the improvement in simulation capability benefited from the CBP is also greater in dry season.

Discussions
Further discussions are carried out in order to investigate impacts of copula selection and correlation coefficient on uncertainty quantification results, taking the wet season as a typical representative.

Impacts of the Copula Selection
Two other CBPs based on Frank and Clayton copulas are also used to quantify the uncertainty in rainfall-runoff simulations for wet season using the identical procedure for comparison purpose. Figure 4 presents the simulative QQ plot, α-index, π-index and CRPS associated with three copulas. It can be seen that compared with the G-H CBP, the Frank CBP has slightly higher reliability but much lower resolution, while the Clayton CBP has lower reliability and resolution simultaneously. The CRPS value in wet season associated with the Frank CBP and Clayton CBP is worsened (increased) by 17.2% and 31.0% compared with the G-H CBP, respectively. This consequence is agreement with the copula fitting results. Consequently, it is highlighted the significance of selecting the appropriate Copula function in the uncertainty quantification of rainfall-runoff simulations using the CBP.

Impacts of the Correlation Coefficient
Ten from the simulative QQ plot that the CBP underestimates simulative uncertainty in wet season for the scenarios where the Kendall tau correlation coefficient is larger than 0.7, otherwise overestimates simulative uncertainty. With the increasing of the Kendall tau correlation coefficient, the α-index decreases firstly, then increases and finally decreases and reaches its maximum at the actual Kendall tau correlation coefficient (0.760). The π-index increases as the Kendall tau correlation coefficient increases. Moreover, The CRPS value firstly decreases, and then increases as the Kendall tau correlation coefficient increases, and reaches its minimum at the actual Kendall tau correlation coefficient. It is demonstrated that effectively capturing and considering the actual dependence between observed and simulated flows is crucial for the CBP to quantify the uncertainty of rainfall-runoff simulations.

Conclusions
Accurate and reliable rainfall-runoff simulation is important for water resources management. This paper presents an uncertainty quantification case study in rainfall-runoff simulations using the copula-based Bayesian processor (CBP) in Danjiangkou Reservoir basin, China. The seasonality of uncertainty in rainfall-runoff modeling is explored, and impacts of copula selection and correlation coefficient on uncertainty quantification results are investigated. The main conclusions were summarized as follows:  1. In this case study from Danjiangkou Reservoir basin, it is found that the observed and simulated flow follows the Log-normal distribution, and the Gumbel-Hougaard is the best copula to capture the actual dependence between observed and simulated flows. The overall performance of the CBP is satisfactory, which provides a useful tool for estimating the uncertainty of rainfall-runoff simulations. 2. The behavior of the CBP is found to be different in wet and dry seasons. The CBP underestimates simulative uncertainty in wet season, while overestimates in dry season. Compared with wet season, the dry season has higher reliability and greater resolution, which illustrates that the CBP captures the actual uncertainty of rainfall-runoff simulations more accurately in dry season. In general, the performance of the CBP in dry season is superior to the wet season and the simulation capability improvement benefited from the CBP is also greater in dry season. 3. The performance the CBP highly depends on the selected Copula function and considered Kendall tau correlation coefficient. As a result, great attention should be paid to selecting the appropriate Copula function and effectively capturing the actual dependence between observed and simulated flows in the CBP-based uncertainty quantification of rainfall-runoff simulations practice.
It is noted that in this paper the best Copula function selected for observed and simulated flows is Gumbel-Hougaard Copula, along with the performance of the CBP in dry season is superior to the wet season is indeed case-specific. This is principally because different hydrological model also has different performance and behavior in rainfallrunoff simulations. In view of the fact that the most appropriate Copula function and seasonal differences in uncertainty assessment performance rely heavily on the hydrological model adopted, developing the CBP based on multi-model for uncertainty quantification and reducing of rainfall-runoff simulations is a valuable research direction. Besides, promoting the application of uncertainty assessment results in rainfall-runoff simulations to risk decision-making of water resources management is promising in the further study.