Study on Urban Rainstorms Design Based on Multivariate Secondary Return Period

With the rapid urbanization, waterlogging losses caused by rainstorm are becoming increasingly severe. In order to reveal the correlations between rainstorm characteristic elements, and make the calculation of rainstorm return period more reasonable and objective, this study established the joint distribution models of rainstorm elements by using copula theory based on the rainfall data in a Chinese megacity, Zhengzhou. Then their combined design values of primary return period (PRP) and secondary return period (SRP) are derived by the maximum probability method and the same frequency method. Finally, the rainstorm pattern was acquired associated with Pilgrim & Cordery method (PC). The results indicate that the calculation of rainstorm return period (RRP) with SRP is more reasonable than PRP. For same RRP, the rainstorm volume (RV) of “Or” return period type is largest, while the “And” return period’s is smallest, and the RVs of Kendall return period and survival Kendall return period are between them. Concerning Kendall return period, the RVs calculated by the maximum probability method and the same frequency method are pretty close, and their relative deviations are from -5.84% to 4.69%. Compared to “Or” return period, the rainstorm patterns of Kendall return period can reduce the magnitude and investment of the stormwater infrastructure. Moreover, the rainfall with designed rainstorm pattern of survival Kendall return period mainly concentrated before the rain peak in contrast with Kendall return period.


Introduction
In recent years, rainstorms are becoming more frequently affected by climate changes and human activities, which have led to serious urban waterlogging disasters and become the focus of public's concern worldwide (Bae and Chang 2019; Bertilsson et al. 2019;Liao et al. 2019;Lohani et al. 2004;Loke et al. 2021;Yao et al. 2021). Zhengzhou city, a superlarge central city of China, also suffers from urban flooding. A typical extreme rainstorm happened in Zhengzhou city was on 7/20/2021, which resulting in 292 deaths and 1.736 million people affected. Thus, the risk assessment of rainstorm event is significant for the urban flood control and disaster reduction.
Rainstorm return period (RRP) was widely used to evaluate the occurrence risk of rainstorm event in previous studies (Zhang et al. 2012;Chin 2017;Singh et al. 2020). RRP refers to the average time interval of a specific rainstorm event. However, traditional calculation of RRP only considers the rainstorm volume, but for other rainstorm characteristic elements as peak rainfall intensity (PR), pre-peak cumulative rainfall, peak appearance time, and rainstorm duration are not involved, thus the correlations of these characteristic elements are ignored. Actually, the multivariable return period are now gradually applied to the hydrology field, such as drought (Chen et al. 2013;Tu et al. 2016; Van de Vyver and Van den Bergh 2018; Zhou et al. 2019), sea wave (Chen et al. 2017) and flood events Li et al. 2014;Requena et al. 2016), but it is still rare for rainstorm. RV and PR are two mainly considerate elements for a rainstorm event, but rainfall after peak (RAP) can also influence the potential inundation process. Therefore, it is crucial to include RV, PR and RAP in rainstorm risk assessment.
"Or" return period (ORP) as well as "And" return period (ARP) are two common used multivariate return period types and also called as primary return period, but they cannot identify the disaster events' hazard area accurately (Salvadori and De Michele 2004;Salvadori et al. 2011;Huang and Chen 2015). To overcome this shortcoming, (Salvadori et al. 2011(Salvadori et al. , 2013 proposed KRP and SKRP in 2011 and 2013 respectively. Compared with ORP and ARP, the theoretical identification of dangerous domain is improved by KRP and SKRP, thus their calculations are more reasonable. In previous research, most concerns are limited to KRP, while SKRP is less, even be ignored, which results in an incomplete study on SRP. In view of this, we first redefined SRP, that is, both KRP and SKRP were included in SRP, and then combined with copula theory, maximum probability method, a new rainstorm pattern was derived, this is also the innovation of this study. This paper analyzed the characteristic elements of rainstorm events firstly, including RV, PR, RAP, and then using Copula method, the ternary joint distribution model of rainstorm elements was constructed. Moreover, combined rainstorm elements with PRP and SRP are displayed to show their differences. Finally, applied with maximum probability method (MP) and same frequency method (SF), a new rainstorm pattern was derived based on designed rainstorm elements. Overall framework of this paper is shown in Fig. 1.

Study Area and Data Sources
Zhengzhou city is located in the south of the North China Plain and the lower reaches of the Yellow River, and it lies between 112°42' E-114°14' E and 34°16' N-34°58' N. There is little precipitation in spring and winter, but more in hot summer season, and its average annual precipitation is about 625 mm. The downtown area of Zhengzhou City is classified into five regions: Huiji District, Zhongyuan District, Jinshui District, Erqi District and Guancheng Ethnic Minority District. Nowadays, concentrated administrative, commercial and residential buildings exist in these regions, once urban flooding occurs, it will produce huge economy losses and public safety will also be threatened.
In this study, the observed rainstorm data of 13 rain gauge stations from 2011 to 2018 are collected from the Zhengzhou Meteorological Bureau. The location of rain gauge stations is shown in Fig. 2.

Copula Theory
It is not easy to establish a joint distribution model with multiple variables if variables are not independent of each other. The Copula function is a successful tool for constructing multi-dimensional joint distribution when the marginal distribution of each variable is known.
Copula theory was firstly proposed by Sklar (Sklar 1959), which can connect the multivariate marginal distributions with a joint distribution. Four commonly used Copula functions in hydrology field are selected to construct the ternary joint distribution of rainstorm elements in this study, including Gumbel Copula, Clayton Copula, Frank Copula and Gaussian Copula. The formulas of Gumbel Copula, Clayton Copula, Frank Copula and Gaussian Copula are shown respectively as follows: where u, v and w are marginal distribution of variables respectively, θ is a coefficient representing the dependence of variates, Φ is the standard normal distribution function, Σ is the correlation coefficient matrix of Gaussian Copula. W is the integral variable matrix.

Secondary Return Period
Supposing a two-dimensional joint distribution (X, Y), its return period is traditionally defined as the following two types (Shiau 2006): 1) if X or Y exceeds the setting thresholds, the calculated return period is called ORP; 2) if X and Y both exceed the setting thresholds, the calculated return period is called ARP. The ORP and ARP were computed based on formula (5) and (6) respectively, and their schematic diagram is shown in Fig. 3.
where u, v and w are the marginal distribution of X, Y, Z, respectively; C is the joint distribution value. As shown in Fig. 3, the curves of " P = P 1 " and " P = P 2 " represent contour lines of two return periods and P means the joint probability. For any point ( u 1 , v 1 ) on the contour line, the corresponding ( x , y ) can be obtained by x = F −1 (u) and y = F −1 (v) . In addition, selecting different points from the same return period contour will result in different risk regions, that is, there are multiple combination events corresponding to the same return period. Moreover, the definitions of ORP and ARP are presented irrationally. For the variables X and Y, the event of "X exceeding the threshold" or "Y exceeding the threshold" is considered as dangerous, larger return period usually means a smaller danger zone, and the dangerous events with larger return period should be involved in the events with smaller return period, but some events that do not satisfy this rule existed in ORP and ARP. Figure 3a is the sketch figure of ORP, "A" is a point on the line P = P 1 and "B" is a point on the line P = P 2 . If P 1 < P 2 , the dangerous region corresponding to point A should be included in the region of point B, whereas, the existing C* area is obviously unreasonable, the same situation happened in Fig. 3b. This indicates that PRP has a theoretical defect in describing the return period of events (Salvadori and De Michele 2004;Huang and Chen 2015).
According to the Kendall distribution function defined by Nelson (Nelsen 2006), three categories of safe, critical and dangerous are classified by Salvadori (Salvadori et al. 2011) through judging the joint distribution probability. In 2013, Salvadori (Salvadori et al. 2013) further divided the two-dimensional space into three parts using curves formed by pair (x, y) satisfying F(x, y) = t . The KRP and SKRP remedied the PRP's deficiency by improving the identification of dangerous zone.F(x) , KRP and SKRP were computed as following: where K C (t) is the Kendall function value with cumulative probability t (Graler et al. 2013), K C (t) is the survival Kendall function value (Salvadori et al. 2013). Since K C (t) and K C (t) cannot be figured out by analytic formula, thus the Monte Carlo simulation method was used to generate rainstorm elements combinations here based on the joint distribution model, then the KRP and SKRP would be calculated with empirical frequency method.

Correlation Analysis
It is found that about 70 percent waterlogging incidents were caused by short-duration rainstorm events, whose lasting time less than 3 h, through python web crawler. According to the statistics, the rainfall events lasting for 1 h accounted for the largest proportion (45 percent), and in which 80 percent present the single-peak rain pattern. Therefore, the rainstorm events with one-peak of 1 h is focused in this paper. Figure 4 shows the correlations of rainstorm element. And the correlation coefficients between RV and PR, RV and RAP, PR and RAP are 0.8623, 0.7446 and 0.5175, respectively, which indicating that there is a positive correlation among them.

Marginal Distribution
Pearson Type III Distribution (P-III), Gaussian Distribution, Gamma Distribution, Generalized Extreme Value Distribution and Weibull Distribution are employed to select the best fitting theoretical marginal distribution. Firstly, Kolmogorov-Smirnov (KS) test was used to proof which can pass the significance test at 0.05 level, and then the fitting results between theoretical and empirical frequency of univariable were calculated. Eventually, P-III distribution was chosen as the optimal type of RV and PR, and Gamma distribution as the best one for RAP. The fitting degree between theoretical and empirical distribution of RV, PR and RAP are 0.9769, 0.9836 and 0.9718, respectively.

Joint Distribution of Rainstorm Elements
The maximum likelihood estimation method (Xu et al. 2008) was applied to calculate the parameters of Gumbel, Frank, Clayton and Gaussian Copula. Furthermore, the square Euclidean distance (d 2 ) ) and OLS value (Kong et al. 2020) of different Copulas are shown in Table 1. The results indicate that the Gumbel Copula performs best for the bivariate joint distributions model of RV-PR and RV-RAP, while for PR-RAP, the Gaussian Copula is the optimal. For the tri-variate joint distribution of RV-PR-RAP, the Gaussian Copula was chosen.

Rainstorm Return Period Calculation
The planning rainstorm return period in the central area of Zhengzhou city is generally within 1 ~ 5 years, and 10 ~ 20 years in the key areas. Therefore, the return periods of 2 years, 5 years, 10 years and 20 years were selected to calculate the recurrence period of rainstorm elements combination here. According to Sect. 2.2.2, the PRP and SRP calculation results of combined rainstorm elements are shown in Table 2.
According to the risk domain represented by each return period type as well as the nondiminishing property of Copula method, for a given return period, ORP should less than KRP and ARP ought greater than SKRP.
As shown in Table 2, for a given return period, ORP is the smallest and even smaller than the given RRP. Besides, the magnitude relationship between the primary return period and secondary return period is ORP < SKRP < KRP < ARP. As mentioned above, ORP's and ARP's identification of risk ranges are both unreasonable, thus PRP will result in lower or higher standards for flood control and drainage engineering, and further bring out waterlogging disasters or finical wasting. For safety reasons, it is rational to adopt KRP rather than ORP in the case of any rainstorm element exceeds the set standard, and SKRP could also be considered as a suitable candidate if all rainstorm elements exceed the set standard.

Design Values of Combined Rainstorm Elements
In three-dimension space, the combined rainstorm elements with same rainstorm return period will form a curve surface. According to (Salvadori et al. 2013;Tu et al. 2018), there exist a combination ( u m ,v m ,w m ) with the maximum joint probability density f (x, y, z) . This means a combination, which most likely happen, could be found out at specific return period level based on the measured rainstorm data. Consequently, this rainstorm elements combination could be applied to the design of rainstorm pattern. The computed method of ( u m ,v m ,w m ) and f (x, y, z) are shown as follows: where u m , v m and w m are the univariate marginal distribution value of rainstorm elements with maximum occurrence probability. In order to reveal the differences among various rainstorm return period (RRP) types as well as the impact of MP and SF methods, the univariate design values of rainstorm elements and their combined design values using MP and SF methods are shown in Table 3.
RV was selected as the major characteristic element and been analyzed here. The results demonstrate that: 1. For the RRP from 2 to 20 years, the RVs of SRP figured by the MP method are smaller than its univariate value, and the RVs of different RRP is sorted as  The above studied results show that, the designed rainstorm values of a single rainstorm element ignored the correlation between rainstorm variables. Besides, the designed values of combined rainstorm elements with PRP may be lower or higher. Therefore, it would be a good choice to use SRP combined with MP method for rainstorm design.

Rainstorm Patterns
Pilgrim & Cordery (PC) method has been widely used to the designing of rainstorm process (Pilgrim and Cordery 1975;Yan et al. 2021). A new rainstorm pattern was introduced on the base of designed rainstorm elements in Sect. 3.3. With 10 min as a time interval, the rainstorm process of 1 h could be divided into six periods. And the implement steps of the new rainstorm pattern are as follows: Step 1: Referring to the determination of the rain peak position with PC method, the rainfall amount in each period of the rainstorm event was sorted firstly, and the period with larger amount has the smaller serial number.
Step 2: Taking the smallest serial number as the rain peak position, and the peak rainfall of SRP calculated by the MP method listed in Table 3 would be placed at this period.
Step 3: Computing out the proportion of rainfall amount in each period, and then their mean proportions with all rainstorm samples could be figured out.
Step 4: Normalizing the proportions in the pre-peak period and post-peak period separately, then according to Sect. 3.3, the designed rainstorm elements values were successfully allocated to each period. Consequently, the new rainstorm pattern was successfully acquired. The rainstorm patterns based on SKPR and MP method are shown as Fig. 5. It can be seen from Fig. 5c that the RV is mainly concentrated in the pre-peak periods, besides, the rain intensity increases rapidly and the RV is accumulated quickly. By contrast, the rainstorm volume of KRP with MP method is mainly focused after the rain peak (Fig. 5b). Moreover, the rainstorm patterns based on the single rainstorm element with its marginal distribution (Fig. 5a) are similar to KRP's. Although the rainstorm patterns of SRP are presented in Fig. 5, the differences between PRP and SRP are needed to be further exploited. In Fig. 6, (i) represents the rainstorm patterns of ORP and KRP, (ii) shows the rainstorm patterns of ARP and SKRP.
The (i) in Fig. 6 denotes that the rainstorm pattern of KRP is nearly involved in the ORP's. This indicates that the waterlogging disaster caused by the later rainstorm pattern is more severe than the former. Once the waterlogging control and drainage engineering is carried out according to ORP type, it is easier to lead to overestimation of rainstorm element value, thus resulting in unnecessary project scale.
The (ii) in Fig. 6 illustrates that the rainstorm volumes estimated by SKRP are mainly concentrated before the rain peak, while the ARPs' are just the reverse at 2, 5 and 10 years return periods level. The situation is quite different from the 20 years', in which the rain peak is more prominent of ARP type.
In addition, Fig. 6 also provides a proof that a better regularity is existed in the rainstorm patterns of PRP and SRP of 2 and 5 years, but for the higher RRP of 10 and 20 years, the regularity is not obviously. For example, for the RRP of 10 years and 20 years, despite the rainstorm pattern of KRP are involved in the ORPs', the excess RV of 10 years return period mainly distributed after the rain peak, while the 20 year's is in reverse. This phenomenon may be explained by the regional rainstorm characteristics and further indicate that rainstorm event with various magnitudes will have their particular rainstorm patterns.
Therefore, compared with ORP, the KRP can decrease the scale of waterlogging control and drainage project by guaranteeing the reasonable return period, thus the project cost will be reduced moderately. The RV of SKRP type is slightly larger than that of ARP, and the RV of SKRP type is mainly concentrated before the peak, while the RV of ARP is in reverse.
Considering the limited rainstorm samples, we selected the rainstorm patterns derived by SRP and MP methods at the return period level of 5 years and 10 years as well the corresponding samples (Fig. 7), and their relative mean deviations are 0.22 and 0.27 separately. Moreover, since the designed combination has not yet occurred, these samples are only close from the RV, whereas there are larger deviations in PR and RAP, so the actual fitting effect could be better.

Conclusions
In this study, the concept of SRP was redefined firstly, and then using Copula theory and MP method, the combination values of rainstorm elements were calculated and the rainstorm patterns were designed successfully.
It is concluded that the extra-threshold risk (shown as Fig. 3) of rainstorm element can be described well by the designed value of combined rainstorm elements calculated by the SRP and MP method, and the designed RV of different return period types is sorted as RV ARP < RV SKRP < RV KRP < RV univariate < RV ORP .
Compared with ORP, the rainstorm patterns derived by KRP and MP method can achieve a proper magnitude of waterlogging control and drainage engineering, thus reduce the unnecessary investment. While compared with ARP, the RV of the designed rainstorm pattern based on SKRP is mainly concentrated before the rain peak, consequently, it would cause different urban waterlogging process and inundation regions. And the comparison between designed rainstorm pattern and samples also indicate the reliability of our research.
However, this study still has some disadvantages to be overcome. For example, the rainstorm series are not long enough, which may lead to less consideration of past rainstorm characteristics. Therefore, if the study is carried out in other areas, more sufficient rainstorm series will have more satisfying results.