Regional Frequency Analysis of Observed Sub-Daily Rainfall Maxima over Eastern China

Extreme climatic events, such as heavy rainfall, floods and droughts, can have severe impacts on society and the economy. Existing approaches to estimating the return periods of extreme events, which furnish a theoretical basis for the planning of hydrological projects, are mainly based on statistical analyses of long-term observational data. Therefore, frequency analysis is of considerable importance for coping with heavy rainfall events and resulting floods (Milly et al., 2008; IPCC, 2014; Jiang et al., 2015; Fischer and Knutti, 2015; Serinaldi and Kilsby, 2015).

Studies have shown that the frequency of heavy rainfall events has clearly increased, while those of light and moderate rainfall have decreased globally in the warming climate (Groisman et al., 2005; Goswami et al., 2006; Moore et al.,2015; Zheng et al., 2016). Similar findings have been reported in numerous studies in China (Liu et al., 2005; Qian et al., 2007; Zhu et al., 2009; Zhang and Zhai, 2011; Jiang et al., 2013; Jiang et al., 2014; Liu et al., 2015). In eastern China, where the East Asian summer monsoon is an important climate feature, atmospheric moisture comes mainly from the western Pacific and Indian oceans. Therefore, the intensity of the western Pacific subtropical high has a crucial effect on the occurrence and spatial distribution of heavy rain and consequent floods (Fischer et al., 2012; Chen et al., 2013; Zou and Ren, 2015). The spatial distribution of heavy rainfall is complex and it is necessary to classify homogeneous regions and select appropriate distributions in each for the frequency analysis of heavy rainfall events. Multivariate techniques, e.g., cluster analysis and rotated EOF analysis, are widely used to classify homogeneous regions of rainfall events (Ramos, 2001; Cowpertwait, 2011).

Based on a probability-weighted moment parameterization method (Greenwood et al., 1979), (Hosking, 1990) and Wallis et al. (1997) introduced the L-moments approach and developed regional frequency analysis (RFA). RFA uses the L-moment of an extreme series to test discordant stations in a basin, and then combines this with a cluster algorithm to divide a large basin into smaller homogenous regions. RFA has been used in many places, including the United Kingdom, United States, Canada, India and China, and has been shown to perform well in regions of data shortage (Hosking and Wallis, 1997; Wallis et al., 2007; Mladjic et al., 2011; Hassan and Ping, 2012; Kumar et al., 2015; Ma et al., 2015). However, most studies have examined homogenous regions at a daily scale within a small region.

The return period of heavy rainfall events is a statistical measurement that typically uses an appropriate parameterization method to fit a candidate distribution, based on annual maximum (AM) or peaks-over-threshold (POT) extreme samples (Coles, 2001). The POT method has a limitation with regard to choosing a proper threshold, so the AM method is used more frequently in related studies (Madsen et al., 1997a, b). There are several types of statistical distributions for frequency analysis of extreme hydroclimatic events. The WMO has recommended the generalized extreme value (GEV) distribution to fit AM series and the generalized Pareto (GPA) distribution for POT series (Klein Tank et al., 2009). Many other distribution families can be candidates for fitting to time series of hydroclimatic extremes. For example, the generalized logistic (GLO) distribution is commonly used in European countries; the Log-Pearson type III (LP-III) distribution is recommended in the United States and Australia; the Gumbel distribution is used in Canada and India; and the lognormal (LN) and the Wakeby distributions have been applied in Japan and Korea, respectively (Park et al., 2001; Öztekin, 2007; El Adlouni et al., 2008; Haddad and Rahman, 2011).

In China, the Pearson-III (P-III) distribution recommended by the Chinese Ministry of Water Resources (1994) has been popular in frequency analysis of floods for engineering design since the 1960s. In addition, there are many other distributions used for frequency analysis of hydroclimatic events. For instance, (Feng et al., 2007) applied the GEV distribution to daily and multi-day (2-, 5- and 10-day) AM rainfall series from China during 1951-2000. They concluded that a 50-yr event of the 1950s had become a 25-yr event in the 1990s in the Yangtze River basin and northwestern China. (Su et al., 2009) used the GEV, GPA, GLO and Wakeby distributions to study AM series of rainfall/drought frequency in the Yangtze basin during 1960-2005, and made future projections using the ECHAM5/MPI-OM model. (Yang et al., 2010) used generalized Gaussian (GG), GLO, GEV and P-III distributions to analyze single and multi-day (3-, 5- and 7-day) AM rainfall series during 1960-2005 in the Pearl River basin. (Fischer et al., 2012) compared the P-III, GEV, GPA and Wakeby distributions in the Pearl River basin, finding that the Wakeby was the best fit distribution in the Pearl River basin, although it had greater uncertainty in parameter estimation.

Many studies have fitted distribution functions to daily data in different river basins. However, disasters are generally caused by hourly rainfall, and there have been limited studies on the frequency of sub-daily rainfall in a large heterogeneous region. In the present study, an RFA approach was used to classify homogeneous regions in eastern China based on AM 1-, 2-, 3-, 6-, 12- and 24-h precipitation series from 1960 to 2014. In the resultant homogeneous regions, return periods of heavy rainfall events were estimated using the optimal distribution of each region and absolute and relative differences of return levels between selected optimal distributions, and the commonly used P-III distribution was then analyzed. The results may be useful to improve capabilities for coping with increased heavy rainfall and consequent floods, and to reduce potential meteorological risks.

In this study, observed hourly rainfall records from 737 meteorological stations across eastern China, spanning the period 1960-2014, were used. The hourly rainfall dataset was examined via homogeneity testing and adjustments, including changes in instrument type, station relocation, and trace biases (Ren et al., 2010). From the original list, 442 stations over eastern China were selected on the condition that less than 5% of data were missing in the rainy season (May-September). The spatial distribution of the selected stations is shown in Fig. 1. The annual maxima of the 1-, 2-, 3-, 6-, 12- and 24-h rainfall series were established using a moving window approach.

Figure 1.  Locations of meteorological stations across eastern China.

2.2.1. Regionalization of homogenous regions by RFA

The RFA based on L-moments was used in this study. This method assumes that the same population L-moment ratios can be found in a homogeneous region for all sites. However, due to the sampling variability, sample L-moment ratios may not be the same (Hosking and Wallis, 1997).

Supposing that there are N stations in a region, each station has an sample of size n, arranged in ascending order, {x₁,x₂,…,x_n}, the L-moments are defined as \begin{equation} \left\{ \begin{array}{l} {l_1=b_0}\\ {l_2=2b_1-b_0}\\ {l_3=6b_2-6b_1+b_0}\\ {l_4=20b_3-30b_2+12b_1-b_0} , \end{array} \right.(1) \end{equation} where \begin{eqnarray} \label{eq1} b_0&=&\dfrac{1}{n}\sum_{i=1}^n {x_i} ,(2) \\ \label{eq2} b_1&=&\dfrac{1}{n}\sum_{i=2}^n\dfrac{i-1}{n-1}x_i ,(3) \\ \label{eq3} b_2&=&\dfrac{1}{n}\sum_{i=3}^n\dfrac{(i-1)(i-2)}{(n-1)(n-2)}x_i ,(4) \\ \label{eq4} b_3&=&\dfrac{1}{n}\sum_{i=4}^n\dfrac{(i-1)(i-2)(i-3)}{(n-1)(n-2)(n-3)}x_i .(5) \end{eqnarray}

The L-moment ratios are defined as \begin{equation} t_2=l_2/l_1,t_3=l_3/l_2,\cdots,t_r=l_r/l_2\quad (r=3,4,\cdots) ,(6) \end{equation} where t₂ is the L coefficient of variation (L _cv), t₃ is L skewness (L _sk) , and t₄ is L kurtosis (L _ku) (Hosking, 1990; Hosking and Wallis, 1997).

A D-test was used to discover the discordant stations in each region. If D_i>3 where a region has more than 15 stations, a station is considered discordant (Hosking and Wallis, 1997). D_i is defined as \begin{equation} \label{eq5} D_i=\dfrac{1}{3}N(u_i-\overline{u})^{\rm T}A^{-1}(u_i-\overline{u}) ,(7) \end{equation} where \begin{eqnarray} u_i&=&(t_it_{i,3}t_{i,4})^{\rm T} ,(8) \\ \label{eq6} \overline{u}&=&\dfrac{1}{N}\sum\limits_{i=1}^N{u_i},(9)\\ \label{eq7} A&=&\sum_{i=1}^N(u_i-\overline{u})(u_i-\overline{u})^{\rm T}. (10) \end{eqnarray} The heterogeneity measures for a region are based on an H-test. Indexes H_1-3 are weighted standard deviations of L _cv, L _sk and L _ku, respectively. In general, H₁ is the most discriminatory of the three indexes. A Monte Carlo simulation can be used to derive H_1-3: \begin{equation} \label{eq8} H=\dfrac{V-\mu_v}{\sigma_v} , (11)\end{equation} where V means the weighted standard deviation of the at-site sample L _cv, which is defined as \begin{eqnarray} \label{eq9} V&=&\left\{\sum_{i=1}^N{n_i(t_i-t_r)^2}\left/\sum_{i=1}^N {n_i}\right.\right\}^{\frac{1}{2}} ,(12)\\ \label{eq10} t_r&=&\sum_{i=1}^N{n_it_i}/\sum_{i=1}^N{n_i} . (13)\end{eqnarray} where, σ_v, μ_v is the mean and standard deviation of V.

If H<1, the region is regarded as homogenous; if H is between 1 and 2, it is possibly heterogeneous; if H>2, it is heterogeneous, and a heterogeneous region can be further divided into smaller subregions by applying a K-means cluster algorithm (Hosking and Wallis, 1997).

2.2.2. Fitting distribution and goodness-of-fit test

Classification of homogeneous regions and parameterization of distributions are both using L-moment based on the statistical characteristics of the sub-daily AM series, and distribution of hourly extreme precipitation within each homogeneous region can be described by an optimal function.

GEV, GPA, GLO, P-III, LP-III, LN and GG distribution functions were used to fit appropriate statistical distributions of extreme values to the sub-daily rainfall time series and estimate the return period in each homogeneous region. These functions are commonly used in studies of daily rainfall extremes. Their probability distribution function f(x) and cumulative distribution function F(x) are shown in Table 1. The parameter estimation of the distributions was based on the L-moments (see details in Hosking, 1990; Hosking and Wallis, 1997).

The Z-test was used to evaluate the goodness-of-fit (GOF) of candidate distributions, as follows: \begin{equation} Z_{\rm Dis}=(\tau_{4,{\rm Dis}}-t_{4, r}+B_4)/\sigma_4 , (14)\end{equation} where, τ_{4, Dis} represents the L _ku of expected distribution, t_4,r represents the L _ku of observed distribution.

Fit a kappa distribution to the regional average L-moment ratios. Simulate a large number, N _Sim, for the mth simulated region, the B₄ representing the bias of t_4,r is: \begin{equation} \label{eq11} B_4=N_{\rm Sim}^{-1}\sum_{m=1}^{N_{\rm Sim}}(t_{4, m}-t_{4, r}) , (15)\end{equation} The standard deviation of t_{4, r} is: \begin{equation} \label{eq12} \sigma_4=\left\{(N_{\rm Sim}^{-1}-1)^{-1}\left[\sum_{m=1}^{N_{\rm Sim}}(t_{4, m}-t_{4, r})^2-N_{\rm Sim}^{-1}B_4^2\right]\right\}^{1/2}. (16)\end{equation}

Figure 2.  Distributions and trends of annual maxima sub-daily rainfall (units: mm) from 1960 to 2014: (a-f) AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall, respectively. Color fill represents the multi-year averaged AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall. Triangles show trends at different levels.

If |Z _Dis| is <1.64, the fitted distribution is acceptable. The distribution with smallest |Z _Dis| is considered the optimal one. The same distribution type was used for calculating the return period of all stations in a given homogenous region.

2.2.3. Mann-Kendall trend test

The Mann-Kendall (M-K) method(Mann, 1945;Kendall, 1970) recommended by the WMO was applied in this study to check for trends in the hydrological series (Du et al., 2014).

To better understand the extreme hourly precipitation frequency over eastern China, here we briefly discuss the long-term climatology of spatiotemporal variation in annual maxima precipitation. Figure 2 illustrates the spatial pattern of long-term (1961-2014) mean AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall, and the nonlinear trends examined by the MK method in eastern China.

The AM 1-24-h rainfall decreased from southeast to northwest. The minimum values of multi-year averaged AM rainfall were found in Inner Mongolia, with the value ranging from 18 mm for 1-h to 52 mm for 24-h rainfall. Maximum values were found in South China, varying between 59 mm and 249 mm for 1-h to 24-h rainfall. Other large rainfall centers were found in the Huaihe, the Sichuan Basin, and the middle/lower Yangtze River, with AM rainfall ranging from 40 mm to 150 mm for AM 1-h durations.

The trends of 1-24-h AM rainfall from 1961 to 2014 showed an upward tendency at 70%-74% of stations over eastern China, but had significant trends at the 95% confidence level at only 15%-20% of stations. The stations where AM sub-daily rainfall had a downward trend were mostly located in northern China, and less than 1% of these stations had significant AM sub-daily rainfall levels at the 95% confidence level. The stations where AM sub-daily rainfall showed obvious changes during 1961-2014 were irregularly distributed in eastern China.

Based on the AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall series during 1960-2014, regional homogeneity tests based on L-moments were used to classify the study area of eastern China into statistically homogeneous regions (Fig. 3). According to the L _cv calculated from sub-daily rainfall extremes, eastern China was divided into six regions: Northeast (NE), Central North (CN), Central (C), Central East (CE), Southeast (SE), and Southwest (SW). The H-test was used in each region to examine statistical homogeneity further. Where the H statistics were >2 for at least two of the six cases (1-, 2-, 3-, 6-, 12- and 24-h rainfall), then further subdivision was performed using a K-means cluster analysis algorithm. As a result, two regions (C and SW) passed the H-test for AM 1-24-h rainfall extremes. The remaining four regions (NE, CE, CN and SE) were subjected to further cluster analysis to obtain smaller homogeneous sub-regions, i.e., NE1 and NE2; CE1, CE2, and CE3; CN1 and CN2; and SE1, SE2, SE3 and SE4. The H statistics for each region are shown in Table 2. It is clear that the H statistics were >1 in some regions for one or two cases, but they were <1 for at least four of the six cases in all 13 regions. Therefore, the heterogeneity test indicates that the 13 regions can reasonably be treated as homogeneous for AM sub-daily rainfall.

Figure 3.  Homogenous regions in eastern China. Color fill represents the L _cv of 6-h rainfall. Solid lines outline regions distinguished by L _cv: (1, 2) Northeast; (3, 4) Central North; (5) Central; (6-8) Central East; (9-12) Southeast; and (10) Southwest. The Northeast, Central North, Central East and Southeast regions are divided into subregions (dashed lines): (1) Northeast 1; (2) Northeast 2; (3) Central North 1; (4) Central North 2; (6) Central East 1; (8) Central East 2; (6) Central East 3; (9) Southeast 1; (8) Southeast 2; (7) Southeast 3; and (12) Southeast 4.

The performance of the distribution functions in each homogenous region was tested by the GOF Z statistics. If a distribution showed good performance for at least four of the six durations (1-, 2-, 3-, 6-, 12-, 24-h rainfall), it was chosen as the optimal regional distribution.

The distributions that passed the Z-test for fitting the single and multi-hour AM rainfall series in each region are shown in Table 3. In the NE1 region, both GEV and GLO distributions passed the GOF Z-test for the 1-h series, but the GEV and GG distributions were appropriate for multi-hour series. Therefore, the GEV was considered the optimal distribution for NE1. For AM rainfall series at all durations, the optimal distribution for NE2 was always the GEV. In CN2, CE1, CE2, SE2, C and SW, the GEV was chosen as the optimal distribution, but it was not the best fit for all durations. In the CN1 region, the GEV distribution best fitted the AM 1-3-h series. However, the GLO distribution was appropriate in five of six cases, so it was selected as the optimal distribution. In SE1, SE3, SE4 and CE3, the GLO was also chosen as the parent distribution because it passed the Z-test for the AM rainfall series for most durations. In summary, the GEV was the optimal distribution for regions NE1, NE2, CN2, CE1, CE2, SE2, C and SW, and the GLO was the optimal distribution for CN1, CE3, SE1, SE3 and SE4.

The 20-, 50- and 100-yr return period events estimated by the optimal distributions of 1-h to 24-h rainfall are presented in Fig. 4. The spatial pattern of return levels for single and multi-hour rainfall series was similar to that of the AM series, with larger values in coastal regions and smaller values in the interior. For a 20-yr event, maximum return levels of AM 1-h to 24-h rainfall were found in the SE4 region, with the value ranging from 80 to 270 mm. The minimum return levels were found in the CN1 and NE1 regions, with return levels ranging from 37 to 104 mm. For the SE1, SE3 and CE1-3 regions, the return levels were about 62-200 mm; for the SE2 and C regions, they were about 60-180 mm; and for the CN2, NE2 and SW regions, they were about 45-128 mm, respectively. The SE1 region had the second largest return period for the AM 1-h to 6-h scale; however, when the time scale was over 12-h, the regional average return level of CE2 and SE3 became larger than that of the SE1 region.

Figure 4.  Spatial distributions of (1) 20-, (2) 50- and (3) 100-yr return levels (units: mm) estimated by the regional optimal distribution: (a-f) represent return levels at 1-, 2-, 3-, 6-, 12- and 24-h AM rainfall, respectively.

Figure 4.  (Continued)

Figure 5.  Absolute differences (units: mm) in the (1) 20-, (2) 50- and (3) 100-yr return levels estimated by the regional optimal distributions and the Pearson III distribution, respectively: (a-f) represent differences in AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall, respectively.

Figure 5.  (Continued)

Figure 6.  Relative differences (units: %) in the (1) 20-, (2) 50- and (3) 100-yr return levels estimated by the regional optimal distributions and the Pearson III distribution, respectively: (a-f) represent the differences in AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall, respectively.

Figure 6.  (Continued)

Figure 7.  Regional absolute differences (units: mm) in the 20- (red), 50- (green) and 100-yr (blue) return levels estimated by the regional optimal distributions and the Pearson III distribution, respectively: (a-f) represent the differences in AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall, respectively; (1-13) represent NE1, NE2, CN1, CN2, CE1, CE2, CE3, C, SE1, SE2, SE3, SE4, SW. The red plus sign means the outliers of the interquartile range.

Figure 8.  The same as Fig. 7, but for regional relative differences (units: %).

The spatial distributions of 50-100-yr return period events showed similar patterns to those from the 20-yr return periods, but with larger magnitudes. From 1-h to 24-h rainfall, maximum return levels were located in SE4, with values ranging from 95 to 333 mm, and from 108 to 390 mm for 50- and 100-yr events, respectively. The minimum return levels were found in CN1 and NE1, with values varying from 45-120 mm and 53-140 mm for 50-, and 100-yr events, respectively. For 50-yr events, the return levels from 1-h to 24-h rainfall in SE1, SE3 and CE1-3 were about 75 mm to 240 mm; for the SE2 and C regions, they were about 70-210 mm; and for the CN2, NE2 and SW regions, they were about 55-140 mm, respectively. For 100-yr events, the return levels from 1-h to 24-h rainfall in SE1, SE3 and CE1-3 were about 85-282 mm; for the SE2 and C regions, they were about 77-238 mm; and for the CN2, NE2 and SW regions, they were about 61-178 mm, respectively.

Both the absolute and relative differences of return levels between the derived optimal distributions and the P-III distribution were compared for 20-, 50- and 100-yr periods. Positive/negative differences mean the return level estimated by the optimal distribution is larger/smaller than that estimated by the P-III distribution. The relative difference is the percentage of the absolute return level differences divided by the return levels estimated by the optimal distribution.

The at-site differences of 20-, 50- and 100-yr return levels for AM 1-24-h rainfall are shown in Figs. 5 and 6, and the regional differences are shown in Figs. 7 and 8. From 1-h to 24-h rainfall, the absolute differences increased with the time scale (1-, 2-, 3-, 6-, 12- and 24-h rainfall), while the relative differences did not. The magnitudes of the 20-yr return levels estimated by the derived optimal distribution were generally smaller than those estimated from the P-III distribution. Over eastern China, the ranges of the return level differences were -3 to 4 mm (-8% to 8%) and -11 to 6 mm (-11% to 4%) for 1-h and 2-h AM rainfall, respectively. The return level differences for 3-h rainfall varied between -8 and 7 mm (-7% to 11%), and that for 6-h rainfall varied between -18 and 6 mm (-11% to 10%). The return level differences increased to -21 to 9 mm (-21% to 9%) and -23 to 11 mm (-13% to 10%) for the 12-h and 24-h AM rainfall.

Return levels for 50- and 100-yr periods of AM 1-h rainfall, estimated by the optimal distributions, were generally larger than those from the P-III distribution. In eastern China, return level differences for AM 1-h rainfall ranged from -4 to 13 mm (-6% to 22%) and -6 to 26 mm (-7% to 43%) for 50- and 100-yr events, respectively, while AM 24-h rainfall ranged from -21 and 60 mm (-7% to 46%) and -10 to 133 mm (-4% to 85%) for 50- and 100-yr events, respectively.

At a regional scale, negative differences for return levels for AM 1-h rainfall were found in the NE2, CN1-2, CE3, C and SE1-4 regions, while positive differences were found in the NE1, CE1-2 and SW regions (Figs. 7 and 8). For AM 2-h rainfall, negative return level differences were found in NE1-2, CN1, CE1-3 and SE1, 3 and 4, and positive differences were found in SE2 and SW. For AM 3-h rainfall, negative differences of return levels appeared in NE2, CN1 and 2, CE1 and 3, C, and SE1-4. For AM 6-h, 12-h and 24-h rainfall, return levels over central and southeast regions (CN1, CE3 and SE1, 3 and 4), estimated via derived optimal distributions, were smaller than those estimated by the P-III distribution.

The spatial pattern of 50- and 100-yr return levels was different from that of the 20-yr return levels. In absolute terms, for AM 1-h to 6-h rainfall, the stations with negative differences of return levels were mostly located in the NE2, CN2, CE1 and 3, C and SE2 regions. For AM 12-h and 24-h rainfall, most of eastern China had positive differences, although the northern CE1 region and part of the SE3 region also had some negative-difference stations.

The RFA of hydroclimatic extreme rainfall can provide a foundation for regional hydrological management, risk assessment and disaster prevention. In this study, regional frequency analyses were conducted on the AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall during 1960-2014, to divide eastern China into homogeneous regions. An optimal regional distribution was selected by a goodness fit Z-test to optimally fit the time series in each homogenous region. Both absolute and percentage differences of 20-, 50- and 100-yr return levels estimated by the optimal distributions and the traditional P-III distribution were analyzed.

Minimum values of multi-year averaged AM rainfall were found in Inner Mongolia, with the value ranging from 18 mm to 52 mm for 1-h to 24-h rainfall, respectively. Maximum values were found in South China, with the value varying between 59 and 249 mm for 1-h to 24-h rainfall, respectively. This is very similar to the results for the 95th percentile threshold value of hourly extreme precipitation given by (Zhang and Zhai, 2011) and (Zheng et al., 2016). During 1961-2014, AM 1-24-h rainfall showed upward trends at 70%-74% stations in eastern China, but these were significant only at 15%-20% stations. Among the stations with a downward trend, less than 1% was significant. The stations with significant changes were irregularly distributed over eastern China. The spatiotemporal variation may be due to rainfall mechanisms. The mechanisms of extreme precipitation over eastern China might be related to the mei-yu front (Fu et al., 2016), tropical cyclones (Fischer et al., 2015) and the mesoscale convective systems triggered by land heating (Chen et al., 2013). The relationship between trends of sub-daily rainfall and rainfall mechanisms will be studied in future work.

Based on L-moments and K-means clustering, the RFA approach divided eastern China into 13 regions: Northeast (NE1 and NE2), Central (C), Central North (CN1 and CN2), Central East (CE1, CE2 and CE3), Southeast (SE1, SE2, SE3 and SE4), and Southwest (SW). Previous studies have divided China into several regions based on daily precipitation and several subjective methods (Jiang et al., 2014; Ma et al., 2015; Li et al., 2016). However, these regions are somehow different to the homogenous regions classified by hourly precipitation data using the RFA method in this study, indicating that the results from daily records cannot be used at hourly scales.

The optimal distributions for the 13 homogenous regions were determined by the goodness of fit Z-test. For the NE1, NE2, CN2, C, CE1, CE2, SE2 and SW regions, the GEV distribution was the most appropriate, whereas the GLO distribution was the optimal distribution for the other regions. The spatial pattern of return levels for single and multi-hour rainfall was large in coastal regions, and small in the interior. Maximum return levels of AM 1-h to 24-h rainfall were found in SE4, with value ranges of 80-270 mm, 95-333 mm and 108-390 mm for 20-, 50- and 100-yr events, respectively. The minimum return levels were in CN1 and NE1, with values ranging from 37-104 mm, 45-120 mm and 53-140 mm for 20-, 50- and 100-yr events, respectively. The spatial pattern of return levels in this study is similar to the studies of (Yao et al., 2009) and (Li et al., 2013), who used hourly data from 1991-2005 and 1982-2012, but with smaller return levels. The return levels of AM 24-h rainfall were different to those of AM daily rainfall, especially in the SE regions (Feng et al., 2007; Su et al., 2009, and Fischer et al., 2012; Jiang et al., 2014; Ma et al., 2015).

The differences in return levels estimated by the regional optimal distributions and the traditional Pearson-III distribution were large, especially for 50- and 100-yr events. For 20-yr events, the absolute differences increased with time scale, with values ranging from -3-4 mm to -23-11 mm (-10%-10%) for 1- and 24-h rainfall, respectively. For 100-yr events, the absolute differences ranged from -6-26 mm (-10%-30%) for 1-h rainfall, and increased to -10-133 mm (-10%-90%) for 24-h rainfall. Considering these large differences in estimated return levels, the RFA method should be used to plan flood management and disaster reduction strategies in the future.

Maximum absolute differences were found in the SE regions, while maximum relative differences were found in the CE regions. The regions with small absolute differences in return levels estimated by the optimal and P-III distributions for AM rainfall may have large relative differences. Therefore, more attention should be given to these regions. The differences in return levels estimated by regional optimal distributions and the P-III distribution did not show similar patterns within some regions. The main reason may be that the P-III distribution was not the optimal distribution for any homogeneous region, and did not fit well for every station in eastern China. The P-III distribution failed to describe AM rainfall at stations that showed large differences compared with the return levels estimated by the regional optimal distribution.

The estimation of return levels is likely to involve many uncertainties. The sources of uncertainty can be classified as assumption uncertainty, input uncertainty and methodology uncertainty. The assumption uncertainty of the return period deduction process is whether the hydrological frequency analysis follows the stationary hypothesis. The RFA used here is based on the stationary hypothesis, but this may be doubtful in an environment under climate change (Milly et al., 2008; Serinaldi and Kilsby, 2015), since global warming has changed the East Asian summer monsoon by modulating the activities of the western Pacific subtropical high and the low-level vorticity (Fu et al., 2016). Although, we did not find a statistical characteristic for the nonlinear trends of AM 1-, 2-, 3-, 6-, 12- and 24-h rainfall. However, the optimal distribution may vary and the return levels may change substantially in the future, as (Mladjic et al., 2011) indicated for Canada and (Su et al., 2009) for the Yangtze River, based on daily extreme precipitation. Different treatments of sample methods, parameterizations, and candidate distributions give rise to methodology uncertainty (Kay et al., 2009; Zhu et al., 2015). Therefore, in subsequent studies we will use global and regional climate models to project future change and study the uncertainties of return-level estimation.