Rainfall erosivity mapping over mainland China based on high density hourly rainfall records

Rainfall erosivity represents the effect of rainfall and runoff on the average rate of soil erosion. Maps of rainfall erosivity are indispensable for soil erosion assessment using the Universal Soil Loss Equation (USLE) and its successors. To 10 improve current erosivity maps based on daily rainfall data for mainland China, hourly rainfall data from 2381 stations for the period 1951-2018 were collected to generate the R factor and the 1-in-10-year EI30 maps (available at https://dx.doi.org/10.12275/bnu.clicia.rainfallerosivity.CN.001; Yue et al., 2020). Rainfall data at 1-min intervals from 62 stations (18 stations) were collected to calculate rainfall erosivities as true values to evaluate the improvement of the new R factor map (1-in-10-year EI30 map) from the current maps. Both the R factor and 1-in-10-year EI30 decreased from the 15 southeastern to the northwestern, ranging from 0 to 25300 MJ mm ha h a for the R factor and 0 to 11246 MJ mm ha h for the 1-in-10-year EI30. New maps indicated current maps existed an underestimation for most of the southeastern areas and an overestimation for most of the middle and western areas. Comparing with the current maps, the R factor map generated in this study improved the accuracy from 19.4% to 15.9% in the mid-western and eastern regions, from 45.2% to 21.6% in the western region, and the 1-in-10-year EI30 map in the mid-western and eastern regions improved the accuracy from 21.7% to 20 13.0%. The improvement of the new R factor map can be mainly contributed to the increase of data resolution from daily data to hourly data, whereas that of new 1-in-10-year EI30 map to the increase of the number of stations from 744 to 2381. The effect of increasing the number of stations to improve the interpolation seems to be not very obvious when the station density was denser than about 10·10 km 1 station.

conservation. The Universal Soil Loss Equation (USLE; Wischmeier andSmith, 1965, 1978) and the Revised USLE (RUSLE; Renard, 1997;USDA-ARS, 2013) have been widely used to estimate soil erosion in at least 109 countries over the past 40 years (Alewell et al., 2019). Rainfall erosivity is one of the factors in the USLE and RUSLE to represent the potential ability of rainfall and runoff to affect soil erosion. In the USLE, erosivity of a rainfall event is identified as the EI value, also denoted as EI30, which is the product of the total storm energy (E) and the maximum 30-min intensity (I30) (Wischmeier, 1959). The 35 erosivity factor (R factor) in the USLE is the average annual total EI values of all erosive events. To recognize interannual rainfall variability, rainfall data of long periods are required (Wischmeier and Smith, 1978). In the original isoerodent map generated by Wischmeier and Smith (1965), stations with rainfall data of at least 22 years were used. To use the USLE, two additional input parameters are required. One is the seasonal distribution of the R factor. To acquire soil erodibility factor (K factor) and cover-management factor (C factor), the seasonal distribution of EI (monthly, Wischmeier and Smith, 1965;or 40 half-month percentage of EI, Renard, 1997, Wischmeier andSmith, 1978) were needed. The other is the 1-in-10-year storm EI value needed to compute the support practice factor (P factor) for the contour farming (Renard et al., 1997).
Calculation of the rainfall erosivity factor requires high-resolution rainfall data to compute E and I30 accurately. In the original study of event rainfall erosivity, the recording-rain-gauge chart was used (Wischmeier and Smith, 1958). Meteorological stations with this type of data are sparsely distributed and the record length is usually quite short. Methods to estimate rainfall 45 erosivity based on more readily available data have been developed widely, such as daily (Angulomartí nez and Beguerí a, 2009; Bagarello and D'Asaro, 1994;Capolongo et al., 2008;Haith and Merrill, 1987;Richardson et al., 1983;Selker et al., 1990;Sheridan et al., 1989;Xie et al., 2016;Yu et al., 1996;Yu and Rosewell, 1996a;Zhang et al., 2002), monthly (Arnoldus, 1977;Ferro et al., 1991;Renard and Freimund, 1994), and annual rainfall (Bonilla and Vidal, 2011;Ferrari et al., 2005;Lee and Heo, 2011;Yu and Rosewell, 1996b). Yin et al. (2015) evaluated a number of empirical models to estimate the R factor using 50 rainfall data of temporal resolutions from daily to average annual, and showed that the most accurate prediction was based on data at the highest temporal resolution.
Once values of the erosivity factor is obtained with site observations, spatial interpolation methods can be used to estimate rainfall erosivity for sites without rainfall data based on surrounding sites to produce the erosivity maps or isoerodent maps.
Breakpoint data and 1-min data are the best datasets for deriving precipitation intensity and estimating rainfall erosivity.
However, 62 stations with 1-min data collected were inadequate for the spatial interpolation of rainfall erosivity over mainland China. More than 2000 stations of hourly data were collected, together with the 62 stations of 1-min data: (a) to develop highprecision maps of the R factor and 1-in-10-year EI30 over the mainland China; (b) to quantify the improvement of the new 70 erosivity maps in higher temporal resolution and station density and better interpolation method compared to current maps.
The study of (Yin et al., 2019) was chosen to represent the latest data had set to estimate the R factor and 1-in-10-year EI30 and related maps. New R factor and 1-in-10-year EI30 maps were produced in this study may improve the estimation of the soil loss in mainland China.  Naipal et al., 20151998-2012 Hourly and sub-hourly 3625 (387 in China) Gaussian Process Regression Panagos et al., 2017Panagos et al., 1980Panagos et al., -1999Panagos et al., 2000Panagos et al., -2017 Daily 30000+ (~800 in China) Thin-plate spline smoothing Liu et al., 2020 *Map of event 1-in-10-year EI30 in China was also generated.

Data
Rainfall data at three temporal resolutions were obtained hourly, daily and 1-min. Hourly data were used to generate maps of 80 R factor and 1-in-10-year EI30. The evaluation of the effect of station density and interpolation methods on generated maps was also based on the erosivity calculated with hourly data. Daily rainfall data were used to adjust the erosivity estimated by hourly data where the record length was short. Data at 1-min interval were used in two ways. One was to calculate the R factor as true values to evaluate the effect of temporal resolution on erosivity estimation along in comparison to the hourly and the daily https://doi.org/10.5194/essd-2020-370 rainfall data. The other is to develop an exponential model for estimating R factor with the mean annual precipitation. The 85 coefficient of the model was used for adjusting the R factor based on hourly data of shorter periods.
Rainfall data at 1-hour intervals from 2381 meteorological stations over mainland China (Fig. 1) were collected and quality

95
The missing data were handled according to the following criteria: (a) a day with more than 4 missing hours was defined as a missing day; (b) a month with more than 6 missing days was defined as a missing month; (c) a year with any missing month in its wet-season was defined as a missing year. The wet-season for stations north of 32°N was from May to September, and https://doi.org/10.5194/essd-2020-370 for those south of 32°N was from April to October. Missing years were removed and missing hours in the remaining effective years were input in two categories: (a) the missing period is followed by a non-zero record, which recorded the accumulated 100 rainfall amount in the missing period based on data notes; (b) the missing period is followed by zero. In the first case, each missing hour and the following non-zero hour were assigned the average value of the non-zero record in these hours. For the second case, the missing hours were input as zero value.
The daily data were obtained for the same 2381 stations over the period of 1951-2014 ( Fig. 1), and precipitation was measured with simple rain gauges. The data were also collected and quality controlled by the National Meteorological Information Center 105 of China Meteorological Administration. Daily data were collected all year round and the number of effective years ranged from 18 to 54 years. Most of the stations (88%) have data of more than 50 years. An effective year of the daily data was defined as there is no more than one missing month in the year, and a missing month was defined as there are more than 6 missing days in the month. The missing records in the effective years were input as zero value.  The R factor map from Panagos et al. (2017) shown in the discussion part of this study was from the global rainfall erosivity dataset published by Joint Research Centre -European Soil Data Centre (ESDAC).

Rainfall erosivity for stations
Hourly rainfall data were first separated into rainfall storms. A continuous period of >= 6 hours of non-precipitation was regarded as the separation of two rainfall storms (Wischmeier and Smith, 1978). Storms with the amount of >= 12 mm were 120 defined as erosive events (Xie et al., 2000), and were used to calculate the rainfall erosivity factors.
Event rainfall erosivity EI30 (MJ mm ha -1 h -1 ) was the product of the total storm energy E (MJ ha -1 ) and the maximum 30-min intensity I30 (mm h -1 ). The maximum 1-hour intensity H60 (mm h -1 ) was obtained in this study: where r =1,2, …, l means that a storm could be divided into l periods or a storm lasted for l hours, er the unit energy (energy per mm of rainfall, MJ ha -1 mm -1 ), Pr the amount (mm), and ir the intensity (mm h -1 ) of the r th hour (USDA-ARS, 2013).
The R factor (Rhour, MJ mm ha -1 h -1 a -1 ) was the mean annual rainfall erosivity and was obtained by multiplying the R factor from hourly data (Rh) and the conversion factor 1.871 : 130 where i=1, 2, …, N means there are N effective years, and j = 1, 2, …, m means there are m erosive rainfall storms in the i th year.
Due to the variability of rainfall erosivity, stations with less than 22 effective years should be excluded (Wischmeier and Smith, 145 1978). However, stations in western China have limited effective years ( Fig. 1). Once stations with less than 22 effective years are removed, the western stations would be too sparse, which would reduce the accuracy of the rainfall erosivity map. To fill the gap of the insufficient years, daily rainfall data observed by simple rainfall gauges (usually have longer periods) of these stations were used.
The calculation of the R factor considered the following two cases: (1) The effective years of hourly data were no less than 150 those of daily data (37% of the stations); (2) The effective years of hourly data were less than those of daily data (63% of the stations). In the first case, the R factor was calculated directly by hourly data with Eq. (1-5). In the second case, the R factor was firstly calculated by hourly data with Eq. (1-5), then adjusted by the mean annual rainfall calculated by daily data (Zhu and Yu, 2015): where Rhadj was the adjusted R factor, Rhour was the estimated R factor using hourly rainfall data, Pd was the mean annual precipitation of longer period (period of the daily data), and Ph was the mean annual precipitation of shorter period (period of the hourly data).
The exponent 1.481 was calibrated with a power function (Eq. 9) of the mean annual precipitation and rainfall erosivity using 1-min and daily rainfall data of 35 stations in China (Fig. 2). The 1-min and daily data share periods of more than 10 years for 160 these 35 stations.
where Rmin was the R factor (MJ mm ha -1 h -1 a -1 ), and Pm was the mean annual precipitation (mm).
R factor based on 1-min rainfall data was calculated using Eq. (10-11): where I30 was the maximum continuous 30-min intensity. Total storm energy (E) was obtained also using methods from RUSLE2 as Eq. (2-3) with a time increment of 1 min.
Calculation of the 1-in-10-year EI30 also considered two cases: (1) The effective years of the hourly data were no less than 22 years (89% of the stations) and the 1-in-10-year EI30 was estimated by hourly data with Eqs. (6-7); (2) The effective years of 170 the hourly data were less than 22 years, but those of the daily data were no less than 22 years (11%), the 1-in-10-year EI30 was estimated by daily data as follows. Firstly, daily rainfall erosivity was obtained by the following equation developed by Xie et al. (2016): where Pdaily was the daily precipitation (≥10mm), parameter α was 0.3937 in the warm season (May to September), and 0.3101 175 in the cold season (October to April). Secondly, the 1-in-10-year daily erosivity was obtained by calibrating the GEV distribution parameters as Eq. (6-7) and the x in the functions was replaced by the annual maximum daily erosivity. Finally, the 1-in-10-year daily erosivity from daily data was multiplied by a conversion factor of 1.17 to correct the 1-in-10-year daily erosivity to approximate the actual event 1-in-10-year EI30 from 1-min data (Yin et al., 2019). The record length was 22 to 29 for 16 (0.7%) stations, 30 to 39 for 44 (1.8%) stations, 40 to 49 for 216 stations (9.1%) stations, more than 50 for 2105 (88.4%) 180 stations when these adjustments were made.

Spatial interpolation and cross validation
Since there is a good correlation between erosivity factor and the mean annual precipitation, the erosivity maps were obtained using the method of Universal Kriging with the annual rainfall as a co-variable. The annual average rainfall was computed using daily rainfall data of the stations and was firstly interpolated using Ordinary Kriging, and the results were used to conduct 185 Universal Kriging. Both the mean annual precipitation and the erosivity factors were interpolated first for each region https://doi.org/10.5194/essd-2020-370  Fig. 1; Li et al., 2014), and then combined to obtain annual precipitation and erosivity maps over China. Buffer areas were used to avoid the discontinuity in the boundary areas following Li et al. (2014).
To evaluate the efficiency of interpolation models, a leave-one-out cross-validation method was applied in each region.
Symmetric mean absolute percentage error (sMAPE) and Nash-Sutcliffe model efficiency coefficient (NSE) were used for the 190 assessment: where n was the number of stations, Fi was the predicted value at the position of the i th station using data from surrounding stations, Ai was the true value at the i th station. 195

The evaluation of the improvement on the accuracy of the erosivity maps
Current erosivity mapping at national scale in mainland China usually uses daily rainfall data from more than 700 stations.
The R factor and 1-in-10-year EI30 maps of Yin et al. (2019) were taken as references to evaluate the improvement in the accuracy of the erosivity maps generated in this study. To compare the accuracy of the erosivity maps of this study and those of Yin et al. (2019), true values from 1-min data using Eq. (10-11) for 62 stations were used to evaluate the improvement of R 200 factor and those from 1-min data for 18 stations (No.1-18; Fig. 2Figure 2) with more than 22 years were used to evaluate the 1-in-10-year EI30. The values extracted from the two erosivity maps for these stations were compared with the true values calculated with 1-min data. Relative error (%) for each station was calculated as follows: where ( ) is the true value calculated from 1-min data, i = 1,2, …, 62 (18); was the value extracted from the erosivity 205 map of this study and that of Yin et al. (2019); j = 1,2,3,4 represents R factor and 1-in-10-year EI30 maps for the two studies; REij was the relative error (%) of the R factor or 1-in-10-year EI30 of the ith station on the jth map. Considering that the absolute values of RE would be high for stations with smaller R factor values, the relative error of the entire map was expressed as the median absolute value of the RE for all the stations.
The erosivity maps in this study used a different procedure as in Yin et al. (2019)  To evaluate the effect of the temporal resolution on the calculated R factor and 1-in-10-year EI30, hourly and daily rainfall data with the same period as the 1-min data at the 62 stations were used. R factors from hourly data were based on Eq. (1-5), those from daily rainfall data were based on Eq. (12), and those from 1-min data were based on Eq. (10-11). The 1-in-10-year EI30 215 were all calculated by calibrating GEV distribution using Eq. (6-7). Erosivity factors from 1-min data can be regarded as the true value. The relative error was computed for evaluating accuracy. To evaluate the effect of station density, hourly data from 774 stations (used in Yin et al. (2019)) and from 2381 stations (used in this study) were used to generate two separate erosivity maps. R factor and 1-in-10-year EI30 values were compared using leave-one-out cross validation method region by region. The sMAPE was calculated for accuracy assessment. 220 To evaluate the effect of interpolation methods, Ordinary Kriging and Universal Kriging with the mean annual rainfall as the co-variable was applied for the R factor and 1-in-10-year EI30 computed using hourly data from 2381 stations. Both interpolation methods were applied to each of four different regions as shown in Fig. 1 and leave-one-out cross validation results were compared. The sMAPE was also calculated to evaluate the accuracy of interpolated values.

Accuracy evaluation on erosivity maps
Taking erosivity maps (R factor and 1-in-10-year EI30) generated by Yin et al. (2019) as references, this study shows a certain improvement in accuracy ( Fig. 3; Table 2). For the R factor, the values in the map of Yin et al. (2019) were underestimated where the R factor was relatively high, and overestimated where the R factor was relatively low. The improvement was particularly noticeable for western China (R < 1000 MJ mm ha -1 h -1 a -1 ) and the southeastern coastal region (R > 10000 MJ 230 mm ha -1 h -1 a -1 ).
Relative errors of the two maps at the 62 stations are shown in Fig. 4

Erosivity maps and changes comparing with the previous study
The R factor generally decreased from the southeastern part to the northwestern part of China (Fig. 5a), ranging from 0 to 250 25300 MJ mm ha -1 h -1 a -1 . The map of 1-in-10-year EI30 shows a similar spatial pattern with that of the R factor (Fig. 5b), ranging from 0 to 11246 MJ mm ha -1 h -1 . Zero R factor value is found at Turpan, Xinjiang Province, where the mean annual rainfall is only 7.8 mm. The maximum of the R factor (more than 20000 MJ mm ha -1 h -1 a -1 ) is located in the southern part of the Guangxi and Guangdong provinces, along the South China Sea, where the mean annual rainfall is more than 2500 mm.
In addition to the overall trend, some local scale characteristics could be identified in the maps. Taking the R factor map as an 255 example, in the western region, the wetter region in northwestern China was located in the west of Dzungaria Basin and along the Tianshan Mountain, which could be captured on the map.

Figure 6: Changes of the R factor(a) and the 1-in-10-year EI30(b) comparing with the previous study
Comparing with the maps from Yin et al. (2019), the new maps can be quite different at some local areas ( Fig. 6a and 6b). The R factor in the new map was higher for most of the southeastern area, and lower for most of the middle and western areas, 265 especially for the southwestern area (Fig. 6a). The change map of 1-in-10-year event EI30 demonstrated similar pattern with that of R factor, whereas with more negative values in some eastern mountainous areas. Figure 7 shows that the R factor estimated from daily data (Eq. 12) is underestimated when the R value is higher than 10000 270 MJ mm ha -1 h -1 a -1 , and slightly overestimated when the value is lower than 2000 MJ mm ha -1 h -1 a -1 . The model using hourly data improved the accuracy by about 11.1% (median value of relative error) compared to that from daily data (Fig. 7). data were used (Fig. 7). Unlike the R factor, 1-in-10-year EI30 was not noticeably improved with an increase in the temporal resolution from daily to hourly data, probably due to the fact that the 1-in-10-year EI30 values estimated using daily data in 275 Yin et al. (2019) had already been multiplied by a conversion factor of 1.17 to correct the 1-in-10-year daily erosivity to approximate the true 1-in-10-year EI30 from 1-min data.

Contribution of the station density
Interpolation for the Western (W) region had the least NSE compared to others, which may be induced by the sparsity of stations (Fig. 1) and the lower spatial correlation of rainfall. The fitted semivariogram for the R factor in W region had a range of 35 km, whereas the ranges for Mid-western (MW), Northeastern (NE) and Southeastern (SE) regions were 288, 261 and 1235 km, respectively. 285 A comparison of cross-validation between the two maps with different station densities shows that the interpolation with denser stations can improve the accuracy by 2.6% ~ 15.4% for the R factor based on the sMAPE index, and by 1.4% ~ 31.8% for 1-https://doi.org/10.5194/essd-2020-370 in-10-year EI30 (Table 3) when the station density increased. The NSE can increase by 0.016 to 0.109 for R factor. For 10-year EI30, the NSE decreased in region W, MW and SE, and increased by 0.038 in region NE.
For the R factor, in region W, the station density doubled (increased from 36.6 to 21.2 10 3 km 2 1 station), and the accuracy 290 improved by 15.4%, whereas the sMAPE of 53.3% was still high with this increase in station density. In region MW, the station density tripled (from 13.9 to 4.8 10 3 km 2 1 station) and the accuracy was improved by 11.1% based on the sMAPE index from 30.7% to 19.6%. In region NE and SE, the station density tripled and quadrupled, respectively, and the accuracy increased about 2.5%. For 1-in-10-year EI30, the improvement for the 1-in-10-year EI30 was even more (by 31.8% in region W and 19.6% in region MW). The improvement was mainly in western regions, and the station density in the eastern China before 295 the increase is enough to describe the spatial variation of the R factor and the 1-in-10-year EI30. It can be inferred that when there were less than about 10· 10 3 km 2 1 station, the increasing of the site density has little impact on the improvement of the interpolation (Fig. 8).

Contribution of the interpolation method
For the R factor, cross-validation of Ordinary Kriging and Universal Kriging with the mean annual rainfall as the co-variable (Table 4) shows that UK improved the interpolation accuracy by 2.3%-9.0% (sMAPE) compared to OK. In the western region, 305 the NSE increased from 0.285(OK) to 0.599(UK). Therefore, it is better to use UK instead of OK when generating the R factor map, especially in western China where station density was low. For 1-in-10-year EI30, UK improved the accuracy by 0.4%-9.7% (sMAPE). In region W, the accuracy improved by 9.7% and the NSE increased from 0.094(OK) to 0.293(UK).

Discussion 310
This study produced quality R factor and 1-in-10-year EI30 maps with hourly data from 2381 stations over mainland China.
The improvement of the R factor map over previously published R factor map can be contributed to the increase in the temporal resolution from daily to hourly data, whereas that of 1-in-10-year EI30 map to the increase of the station density in comparison with those of Yin et al. (2019). There are mainly two reasons for this. First, 1-in-10-year event EI30 values estimated from the daily data had already been adjusted to those from the 1-min data by multiplying a conversion factor of 1.17 (Yin et al., 2019) , 315 which resulted in no obvious improvement from the daily data to hourly data. Second, the 1-in-10-year event EI30 associated with extreme rainfall event intrinsically has a high spatial variability in comparison to the annual average rainfall erosivity as shown in Table 3. The accuracy of spatially interpolated rainfall erosivity was more sensitive to the station density when the station density is low. Hence the improvement of the map of the 1-in-10-year EI30 was mainly contributed to the increase of the station density, especially for the western and the mid-western regions with sparse station density. 320 Panagos et al. (2017) developed a Global Rainfall Erosivity Database with hourly and sub-hourly rainfall data from 3625 stations over 63 countries, for water erosion assessment for many regions of the world especially where observational data were limited. In their study, rainfall data at 60-min interval from 387 stations across China were used. Figure 9 shows that the R factor for China extracted from Panagos et al. (2017) is systematically underestimated by about 30% for most areas in China, whereas overestimated in the Tibetan Plateau (cf. Fig. 5a). The reason for the underestimation may be that the R factor 325 calculated from 60-min interval data applied a conversion factor (CF30) that was developed from the values estimated by 60min data to those by 30-min data in Panagos et al. (2015), rather than a factor to those by breakpoint data (CFbp) or 1-min data (CF1), which were used in USLE (Wischmeier andSmith, 1965, 1978), RUSLE (Renard, 1997) and this study. Previous research have showed the difference between CF30 and CFbp (CF1) can result in an underestimation of R factor by about 20% (Auerswald et al., 2015;. Table 5 shows that the relative error of the map from Panagos et al. (2017) could 330 reduce by about 6.2% after multiplying by a conversion factor of 1.253, which was calibrated by  for converting the R factor from 30-min data to 1-min data. The adjusted map still generally underestimated. The reason may be that the equation for estimating the storm energy (E) used in Panagos et al. (2017) was from RUSLE (Renard, 1997), which have been reported an underestimation of the storm energy up to 10% in previous studies (McGregor et al., 1995;Yin et al., 2017). Because of this, the equation for estimating the storm energy (E) in RUSLE (Renard, 1997) was then modified in 335 RUSLE2 (USDA-ARS, 2013), which was adopted in this study.
The R factor in the Tibetan Plateau varies from 0 to 12326 MJ mm ha -1 h -1 a -1 in Panagos et al. (2017), and from 4.6 to 4441.9 MJ mm ha -1 h -1 a -1 in this study. The former was derived from a Gaussian Process Regression (GPR) model and a number of https://doi.org/10.5194/essd-2020-370 monthly climate variables from the WorldClim database, such as the mean monthly precipitation, mean minimum, average and maximum monthly temperature. The GPR model was calibrated using the site-specific R factor values and these climate 340 variables, which may not applicable for sites at high altitude, as none of the sites was located in the Tibetan Plateau region.
The GPR model might be the main reason for the overestimation of the R factor in the Tibetan Plateau where the R factor was expected to be underestimated just like any other regions.

Conclusions
This study generated the R factor and 1-in-10-year EI30 maps using hourly rainfall data for the period from 1951 to 2018 from 2381 stations over mainland China. The improvement in the accuracy of these erosivity maps was evaluated against the current 355 https://doi.org/10.5194/essd-2020-370  Yin et al. (2019) were taken as references) in terms of temporal resolution of the rainfall data, the station density, and the interpolation method. The conclusions were drawn as follows: (1) Comparing with the current maps for the 62 reference sites, the R factor map generated in this study improved the accuracy from 19.4% to 15.9% in the mid-western and eastern regions, 45.2% to 21.6% in the western region, and the 1-in-10-year EI30 map improved the accuracy from 21.7% to 13.0% in the mid-western and eastern regions. 360 (2) The R factor and the 1-in-10-year EI30 increased from the northwestern to the southeastern China. The R factor was from 0 to 25300 MJ mm ha -1 h -1 a -1 , and the 1-in-10-year EI30 was from 0 to 11246 MJ mm ha -1 h -1 . Comparing with the current maps, the R factor and 1-in-10-year event EI30 in the new maps were higher for most of the southeastern area, and lower for most of the middle and western areas.
(3) The improvement of the R factor map can be mainly contributed to the increase of the temporal resolution from daily to 365 hourly, whereas that of 1-in-10-year EI30 map to the increase of station density. The increased station density mainly improved the accuracy in the western regions for both the R factor and 1-in-10-year EI30. The contribution of increasing the station density to improve the interpolation was limited when the station density was denser than about 10· 10 3 km 2 1 station. As for the interpolation method, Universal Kriging with the mean annual rainfall as the co-variable performed better than Ordinary Kriging for all regions, especially for the western regions. 370