To verify the superiority of the proposed model in robustness and assessment accuracy, comparative experiments are designed from two aspects: simulating the robustness and fitting accuracy of univariate marginal distributions, and simulating the stability and fitting accuracy of bivariate joint probabilities. Comparative experiments are conducted using Matlab programming. Finally, the experimental results are compared and analyzed to provide evidence and insights.
3.4.1. Comparison Experiments on the Accuracy and Robustness of Marginal Distribution Simulations
To validate the superiority of KDE in simulating the univariate probability distribution of monthly supply–demand sequences, four commonly used distribution functions in the hydrological field (gamma [
52], normal [
48], logistic [
49], Pearson3 [
57]) were selected for comparison experiments with the introduced KDE in this paper. The Kolmogorov–Smirnov test (KS test) and root mean square error (RMSE) were utilized to evaluate the robustness and assessment accuracy of the comparison methods from both confidence probability and fitting accuracy perspectives.
The Kolmogorov–Smirnov test [
71,
72] is a non-parametric goodness-of-fit test. It determines whether to reject the null hypothesis by comparing the maximum absolute deviation (D-statistic) between the empirical cumulative distribution function (ECDF) of the sample and the theoretical cumulative distribution function (CDF) to the critical value at a specified significance level (0.05). Suppose the sample set
has an empirical cumulative distribution function
, and
X follows a theoretical distribution
. The formula to compute the maximum absolute deviation
between the two is as follows:
where
represents the supremum, which is the lowest upper bound of the distances, i.e., the maximum value among all possible absolute differences. If
follows the theoretical distribution
, then as
tends to infinity,
almost surely converges to 0.
The root mean square error (RMSE) [
73] is commonly used to assess the performance of prediction models or the goodness-of-fit of data. In this paper, RMSE is utilized to quantify the difference between the empirical distribution and the fitted theoretical distribution as a measure of goodness-of-fit evaluation [
74]. RMSE can be expressed as follows:
In the above equation, represents the empirical cumulative probability distribution of the observed data,
represents the fitted theoretical probability distribution, and is the number of samples in the dataset.
We simulated the probability distributions of monthly and annual water supply and demand sequences for 13 assessment units (sub-regions) in the study area using gamma, normal, logistic, Pearson3, and KDE, resulting in a total of 169 sequences for both supply and demand. Then, based on the fitted models and empirical cumulative probabilities, we used Equations (12) and (13) to calculate the KS p-value and RMSE.
- ①
Comparative analysis of the Kolmogorov–Smirnov test results
The results of the KS
p-value calculations for the demand and supply sequences are plotted as curves in
Figure 6. The
p-value represents the probability of observing the current sample data or more extreme data under the assumption that the null hypothesis is true. A higher
p-value indicates a higher probability of observing the data, indicating a better fit. From
Figure 6a, it can be observed that out of the 169 demand sequences, only the 164th sequence has a KS test result of 0.038 (normal distribution fit), which is below the significance confidence level of 0.05. The results of the KS test for the other sequences are above the confidence level, indicating that extreme situations in water demand are relatively rare, and most tend towards a normal distribution. From the KS
p-value calculation results for the supply sequences (
Figure 6b), it can be seen that, except for KDE, the KS tests for the other four distributions show instances where the
p-value is below the 0.05 confidence level. This indicates that the distribution of water supply is complex and diverse, and some sequences do not belong to any of the four known distributions; only KDE can correctly fit their probability distribution. Therefore, this indicates that KDE can perform well when the fitting results of other methods are poor. It can adapt to all extreme distribution shapes, demonstrating robustness and resilience unmatched by other methods.
The statistical results of the KS
p-values, including the mean, variance, and rejection rate of the null hypothesis, are summarized in
Table 6. The optimal values are highlighted in bold in the table. The mean is used to measure the overall fitting accuracy of the method. A higher mean of the KS
p-value indicates higher overall fitting accuracy. From the mean statistics in
Table 6, it can be observed that KDE has the highest overall fitting accuracy for both the demand and supply sequences. The variance is used to evaluate the volatility of the KS
p-value results, which reflects the stability of the method. A smaller variance indicates less volatility and greater stability of the method. From the variance statistics in
Table 6, it can be seen that KDE exhibits the most stable fitting of the distribution shapes for both the demand and supply sequences. The rejection rate of the null hypothesis indicates the proportion of
p-values lower than the significance level of 0.05, serving to characterize the fitting capability of the method when facing various distribution shapes. As shown in
Table 6, the rejection rate of the null hypothesis for KDE is 0, indicating that KDE can adapt to the distribution shapes of all water supply and demand sequences, demonstrating the strongest modeling capability.
- ②
Comparative analysis of the RMSE evaluation results
The RMSE calculation results for the empirical and theoretical probability distributions of the water demand and supply sequences are shown in
Figure 7. RMSE is a measure of the deviation between empirical and theoretical distributions, where smaller RMSE values indicate better fitting accuracy. The mean, variance, and range of RMSE values are then calculated in
Table 7, with the optimal values highlighted in bold. The mean RMSE represents the average deviation of the fitting, with smaller values indicating higher overall fitting accuracy. The variance and range of RMSE values assess the stability and robustness of the fitting, where smaller values indicate lower fluctuation and a narrower range. From
Figure 7a and
Table 7, it can be observed that all five methods exhibit relatively high fitting accuracy for the water demand sequences, with fitting deviations below 0.1. Pearson3 and KDE have the smallest mean deviations, but KDE also demonstrates the lowest variance and range, indicating that among the five methods, KDE not only offers high fitting accuracy but also superior stability and robustness. Regarding the RMSE evaluation results for the water supply sequences in
Figure 7b and
Table 7, it is evident that the fitting deviations vary significantly among the five methods. Except for KDE and logistic, the other three methods show considerable fluctuation, especially gamma, which has the largest mean deviation, variance, and range. This suggests that the gamma distribution cannot adequately accommodate all distribution patterns of the water supply sequences. In contrast, KDE exhibits the smallest mean deviation, variance, and range, indicating that KDE not only offers the highest fitting accuracy but also the best stability and robustness, capable of accommodating all extreme distribution shapes of the water supply sequences.
- ③
Comparison and analysis of the fitting results
To visually compare the fitting effects of the five methods more intuitively, empirical probability distribution graphs and fitted theoretical distribution graphs were further plotted. Five representative sets (
Figure 8) were selected from all the results for comparative analysis. From the fitted curves, it can be observed that for sequence data tending towards a normal distribution, such as the water demand and supply sequences of **ushan, all methods achieve relatively good fitting effects. For sequences with significant fluctuations in demand and supply, such as Luohe and Yanhe, except for poor fitting by the gamma distribution, the other four methods can still adapt. However, for supply sequences with special distribution shapes like Shifu and Lishan-2-month, only KDE can accurately fit their probability distributions.
Overall, it is evident that KDE demonstrates the optimal fitting capability, showing robustness and resilience for random variables of water supply and demand with diverse and unknown distribution shapes. It can adapt to extreme distribution shapes and is highly suitable as a marginal distribution simulation function for water supply and demand sequences.
3.4.2. Joint Probability Simulation Accuracy and Robustness Comparative Experiment
The accuracy of marginal probability fitting for water supply and demand univariate variables, along with the suitability of copula functions, jointly determine the accuracy of the joint probability simulation. The higher the accuracy of the joint probability simulation, the more accurate the assessment of the water shortage risk. To further verify the superiority of KDE in joint probability simulation, the marginal probability distributions of water supply and demand obtained by the previous five methods are, respectively, input into the five copula functions in Equation (6) to calculate their respective joint probability distributions. That is, the joint probabilities of water supply and demand are calculated for each combination of the five univariate probability simulation methods and the five copula functions. Then, the mean of the Spearman correlation coefficient and Kendall rank correlation coefficient for the bivariate copula under the given parameters, which are shown in
Table 5, is used to select the optimal copula for each method as the final simulated linkage function. Finally, the accuracy of the model is evaluated by calculating the squared Euclidean distance (SED) between the empirical joint probability distribution and the theoretical joint probability distribution.
The Spearman and Kendall rank correlation coefficients are non-parametric methods for measuring the strength and direction of the relationship between two variables based on the ranks of data objects [
65]. They are particularly suitable for situations where the data do not follow a bivariate normal distribution or the measurement scale is not continuous and quantitative, and they are not influenced by outliers. The difference lies in that the Kendall rank correlation coefficient (
) is based on the concordant and discordant pairs of two sample datasets, as shown in Equation (14), while the Spearman rank correlation coefficient (
) is based on rank differences, as shown in Equation (15).
In the equations,
represents the number of samples,
represents the total number of pairwise combinations of samples,
and
, respectively, represent the number of concordant and discordant pairs;
is used to handle tied ranks, where
and
represent the number of tied ranks in datasets
and
, respectively. It is important to note that tied ranks occurring simultaneously in both
and
are not counted in
and
.
In the equation, represents the number of samples, and is the absolute difference between the ranks of the original observed data and .
The formula for calculating the mean of the Spearman coefficient (
) and Kendall rank coefficient (
) (abbreviated as
MSK) is as follows:
The formula for calculating the squared Euclidean distance between the empirical joint probability distribution and the theoretical joint probability distribution is as follows:
The symbol represents the empirical joint probability distribution, represents the theoretical joint probability distribution fitted by different copula functions, and is the sample size.
- ①
Comparative analysis of MSK
We use Matlab’s Copulastat function to calculate the
and
of the bivariate copula under the given parameters, to measure the suitability of the input data with the copula function. Higher means of
and
(
MSK) indicate that the copula function is more capable of describing the dependence between input data, making it more suitable for simulating their joint probability distribution. Therefore, the copula function corresponding to the maximum
MSK is chosen as the final connection function to compute the joint probability. The calculated results of
MSK for the optimal bivariate copula, using gamma, normal, logistic, Pearson3, and KDE as marginal distribution functions, are shown in
Figure 9. The mean and variance of
MSK are further calculated and presented in
Table 8, with the optimal values highlighted in bold. From
Figure 9, it can be observed that the
MSK of the bivariate copula with KDE as the marginal distribution function is higher than that of the other four methods. Particularly for some extreme supply and demand distribution shapes, KDE demonstrates significant superiority and robustness. From the statistical results in
Table 8, it can be observed that the bivariate copula with KDE as the marginal distribution exhibits the optimal ability to fit joint probability distributions, demonstrating strong stability and robustness, particularly for extreme bivariate distribution shapes.
- ②
Comparative analysis of SED
First, the empirical joint probability distribution values of the supply and demand marginal distributions were calculated, and then, the optimal bivariate copula joint probability distribution values were calculated using five methods, gamma, normal, logistic, Pearson3, and KDE, as the marginal distribution functions. Subsequently, these values were input into Equation (17) to compute the squared Euclidean distance (SED) between the empirical joint probability (ECDF) and theoretical joint probability (CDF). The calculation results are shown in
Figure 10, while the mean and variance statistics of SED are presented in
Table 9, with the optimal values highlighted in bold. From
Figure 10 and
Table 9, it can be observed that the optimal bivariate copula model based on KDE exhibits minimal fitting bias and low fluctuation for all supply and demand sequences. The SED’s mean and variance are both minimal, indicating high simulation accuracy, robustness, and stability of the optimal bivariate copula model based on KDE.