A Universal Multifractal Approach to Assessment of Spatiotemporal Extreme Precipitation over the Loess Plateau of China

Extreme precipitation (EP) is a major external agent driving various natural hazards in the Loess Plateau (LP), China. Yet, the characteristics of spatiotemporal EP responsible for such hazardous situations remain poorly understood. We integrate universal multifractals with a segmentation algorithm to characterize a physically meaningful threshold for EP (EPT). Using 15 daily data from 1961 to 2015, we investigate the spatiotemporal variation of EP over the LP. Our results indicate that, with precipitation increasing, EPTs range from 17.3 to 50.3 mm·d while the mean annual EP increases from 35 to 138 mm from northwest to southeast LP. Further, EP frequency (EPF) has historically spatially varied from 54–116 d, with the highest EPF occurring in the mid-southern and southeastern LP where precipitation is much more abundant. However, EP intensities tend to be strongest in the central LP where precipitation also tends to be scarce, and get progressively weaker as we move towards 20 the margins (similarly with EP severity). An examination of atmosphere circulation patterns indicates that the central LP is the boundary where tropical cyclones reach furthest into inland China, resulting in the highest EP intensities and EP severities being in this area. Under the control of the East Asian monsoon, precipitation from June to September accounts for 72% of the total amount, while 91% of the total EP events are concentrated in June to August. Further, EP events occur, on average, 11 days earlier than the wettest part of the season. These phenomena are responsible for the most serious natural hazards in the 25 LP, especially in the Central region. Spatiotemporally, 91.4% of the LP has experienced a downward trend of precipitation, while 62.1% of the area has experienced upward trends of the EP indices, indicating the potential risk of more serious hazardous situations. The universal multifractal approach considers the physical processes and probability distribution of precipitation, thereby providing a formal framework for spatiotemporal EP assessment at regional scale.


Introduction 30
Extreme precipitation (EP) is the dominant external agent driving floods, erosion, and debris flow etc., with adverse impacts on human life, the social economy, the natural environment, and ecosystems (Min et al., 2011;Pecl et al., 2017;Walther et al., 2002). These impacts are especially severe in arid and semiarid areas because of the sparsity of vegetation and fragility of the eco-environment (Bao et al., 2017;Huang et al., 2016). In recent decades, worldwide climate change has given rise to spatially heterogeneous changes in EP regime (Donat et al., 2016;Manola et al., 2018;Zheng et al., 2015). Such uneven changes in EP 35 have the potential to aggravate adverse impacts on human life and eco-environment and, consequently, EP has recently received increased attention (Eekhout et al., 2018;Li et al., 2017;Wang et al., 2017). In this regard, a fundamental need is to evaluate the regional spatiotemporal variation of EP, providing important information that is crucial to natural resource management and sustainable social development.
The Loess Plateau (LP), which is located in the middle reaches of the Yellow River, is a typical arid and semiarid region, 40 characterized by serious EP-induced natural hazards, including soil erosion and consequent hyper-concentrated flooding and occasional landslides, etc. (Cai, 2001). In the LP area, EP-induced soil erosion generates some of the highest sediment yields observed on earth, ranging from 3 × 10 4 to 4 × 10 4 t·km -2 ·yr -1 . For example, sediment delivered to the Yellow River in recent decades was estimated to be 16 × 10 8 t·yr -1 (Ran et al., 2000;Tang, 2004), and sediment deposition resulted in the river bed of the lower Yellow River aggrading by 8-10 m above the surrounding floodplain (Shi and Shao, 2000). As it flows on the 45 aggraded thalweg, the extreme-precipitation-driven hyper-concentrated floodwaters cause, in turn, the lower Yellow River to burst its channel. Over the past 2500 years, this has cause flooding 1,593 times, and changes to the course of the river channel 26 times, leading to unimaginable death and devastation (Ren, 2006;Tang, 2004). To control such EP induced natural hazards, ecological restoration projects have been implemented over the LP. For example, the "Grain for Green" project (the largest investment project in China) was implemented to control natural hazards such as soil erosion and flooding, and has cost more 50 than 75 billion dollars over the past 20 years.
Accordingly, a better understanding of spatiotemporal EP changes in this area is of considerable interest for various fields, such as risk estimation, land management, flood control, and infrastructure planning (Feng et al., 2016;Wang et al., 2015).
Considerable past work has been devoted to investigating the spatiotemporal variation of total precipitation and precipitation extremes in this region, with consensus obtained for precipitation "amount" (Li et al., 2010a;Li et al., 2010b;Miao et al., 2016;55 Wan et al., 2014;Xin et al., 2009). However, the spatial distribution of EP in the LP is still poorly understood, with considerable disagreement regarding EP and the inability to account for the spatial distribution of EP-induced hazards such as soil erosion (Li et al., 2010a;Li et al., 2010b;Miao et al., 2016;Wan et al., 2014;Xin et al., 2009). It is therefore important to account for the spatiotemporal role of EP in natural hazards, to facilitate better catchment management in regards to issues such as freshwater shortage (Feng et al., 2012). 60 3 To understand spatiotemporal variations in EP, scientists are often required to collect more detailed data, including maximum depth, duration, and area observations (Duliè re et al., 2011;Herold et al., 2017;Miao et al., 2016). Despite sophisticated methodologies, such efforts rely on data from various sources, which are typically absent in the long-term historical observational records, especially over large areas. Therefore, any investigation of spatiotemporal variation in EP must make use of the information in available historical data that was observed at fixed time intervals (e.g. daily). In EP 65 assessment using historical daily data, it is a crucial step to identify EP events by extreme precipitation threshold (EPT) determination. However, EP events tend to be relatively rare, unpredictable, and often of short duration (Liu et al., 2013), and this uncertainty, combined with varying geographical and meteorological conditions, increases the complexity of EP assessment.
In general, existing methods for EPT determination can be grouped into two categories: nonparametric and parametric. Non-70 parametric methods use fixed critical values or percentiles to define the thresholds for extreme events. Because the corresponding classification of EPTs varies from region to region (e.g., a 50 mm daily precipitation event is considered normal in South China but would be an EP event in the LP), the application of non-parametric methods can require considerable subjectivity (Liu et al., 2013), and significantly affect the results of the analysis. For instance, using an absolute value of 50 mm/d, Xin et al. (2009) reported a spatiotemporal decreasing zone of EP in mid-eastern LP while, using the 95% percentile to 75 determine EP, Li et al. (2010a) and Li et al. (2010b) found an increasing trend of EP frequency in the southeastern LP. These reports did not, however, explain the rational for spatial variation of EP and its impacts on the most serious soil erosion in the central LP.
Parametric statistical methods based on empirical distributions have recently become popular. A variety of special distribution functions and parameter estimation techniques have been proposed to characterize observed EP (Anagnostopoulou 80 and Tolika, 2012;Beguerí a et al., 2009;Deidda and Puliga, 2006;Dong et al., 2011;Li et al., 2005;Pfahl et al., 2017). A recent focus that has emerged is to obtain better physical understanding of EPTs, and thereby to assess regional variations in EP. For example, Liu et al. (2013)  Universal multifractals were conceived to study the multiplying cascades governing dynamics of various geophysical fields (Lovejoy and Schertzer, 2013;Schertzer and Lovejoy, 1987). For precipitation, a scaling break separating the meteorological and climatological regimes varies from several days to 1 month, with an average about two weeks (Tessier et al., 1996;Tessier et al., 1994Tessier et al., , 1993. The meteorological scaling interval indicates that, from the multifractal perspective, data collected at time 90 intervals of one day and those at intervals finer than minutes can equivalently characterize the physical processes associated with precipitation (Pandey et al., 1998;Tessier et al., 1996), indicating that EP events can, in principle, be characterized by the 4 study of daily data observed at gauging stations. Of course, it is vitally important that the universal multifractal characterizes how extremes occur in a natural manner (Lovejoy and Schertzer, 2007;Tessier et al., 1996).
In this study, we use the universal multifractal technique to obtain a physically meaningful characterization of EPT. Our 95 objectives are to (1) apply the universal multifractal approach to determine a unique set of EPTs for the LP area, (2) investigate how spatial variations in EP are responsible for the severe nature of soil erosion, and (3) assess the spatiotemporal variation of EP over the LP during the period 1961-2015.

The relationship between precipitation extremes and Multifractals 100
The approach outlined below was used to identify EP events at the observation time scale. In the method of universal multifractals, two equivalent routes can be followed to investigate time series scaling: the probability distribution and the statistical moments. A fundamental property of multifractal fields related to the probability distribution is given by the equation (Lovejoy and Schertzer, 2013;Schertzer and Lovejoy, 1987): where λ represents the resolution of the measure (i.e., the ratio of the external scale L to the measurement scale l; λ = L/l), φλ is the intensity of the field (in this case accumulated precipitation) measured at resolution λ, γ is the order of singularity (maximum precipitation) corresponding to φλ, and the codimension function c(γ) characterizes the sparseness of the γ-order singularities (this function is a basic multifractal probability relation for cascades). Accordingly, the statistical moments are given by 110 where K(q) is the multiple scaling exponent function for moments; q is the order of the statistical moment and   denotes the average of the field (averaged precipitation) at scale ratio λ. The two equivalent routes are related via a Legendre transform (Parisi and Frisch, 1985). The universal K(q) functions and the codimension function c(γ) are expressed as whereα is the multifractal index, which describes how rapidly the fractal dimensions vary as we leave the mean. For time series in this paper, 0 <  <1, and ′ is the auxiliary variable defined by 1/α′ + 1/α = 1 (Lovejoy and Schertzer, 2013). The term C1, the codimension of the mean of the process, varies on 0 ≤ C1 ≤ D (D is space dimension; D = 1 for time series) and quantifies the sparseness of the mean. In this paper, the parameters  and C1 of the multifractal model were estimated using the double 120 trace moment technique (Lavallé e et al., 1993).
As noted by Gagnon et al. (2006) and Lovejoy and Schertzer (2007), the parameters C1 and α characterize the mean of the field, whereas the extremes are expressed by the singularity, γ, and the codimension function c(γ) (Hubert et al., 1993). For 0 ≤ α < 1 and considering a time series of infinite length, a finite maximum order of singularity, γ0, can be determined as However, in general, any time series of finite length will almost surely miss the presence of rare extremes in the field. In this case, the observed singularities will be bounded by an effective maximum singularity, γs: where D is the embedding space dimension (D=1 for the time series). The total dimension of this problem is actually (D+Ds), where Ds=logNs/logλ is the sampling dimension, and Ns is the number of independent time series at each location. The 130 parameter γS links the physical processes that generate precipitation events to the conceptual model of multiplicative cascades, and allows the extremes to be cast in a probabilistic framework, γS > 0. Thus, it was used to infer the extreme events of precipitation fields over well-defined ranges of scale.
For parameter estimation, the parameters α and C1 of the multifractal model are estimated by the double trace moment (DTM) technique (Lavallé e et al., 1993). The q,η DTM at resolution λ and Ʌ is defined as 135 where the sum is obtained over all the disjoint D dimensional balls Bλ,i (with intervals of length τ = T/λ) that are required to cover the time series. K(q,η) is the double trace scaling exponent, and K(q,1) = K(q) is the scaling exponent. This relation can be expressed as where the notation indicates that the multifractal φ at a (finest) resolution Ʌ is first raised to the power η then degraded to resolution λ, and the qth power of the result is averaged over the available data. The scaling exponent K(q,η) is related to K(q,1) ≡ K(q) and given by In the case of universal multifractals, by plugging Eq. (3) into Eq. (5), K(q, η) has a particularly simple dependence on η: 145 where α can be estimated from a simple plot of K(q,η) versus log (η) for fixed q. Thus, based on DTM technique, all parameters can be estimated.

Determination of the extreme precipitation threshold
The approach outlined below was used to estimate the EPTs and EP events for each station, using the singularity parameter, 150 γs. Given that the multifractal parameters α and γs naturally characterize extremes, both of these parameters will change if we gradually remove extreme values from the data set. The singularity of the precipitation data series will be completely changed, and the two parameters will significantly change if all of the extreme values are deleted. To obtain a physically meaningful value for the EPT, we attempted to estimate the multifractal parameter series by applying the universal multifractal approach to our precipitation series and successively eliminating maximum values of precipitation. However, as shown by Figures 1a  155 and 1c, the degree of convergence to the original value is not unique, with the values fluctuating slightly around the original α and γs. Accordingly, the variance series of index series αj and γs,j were computed to eliminate the fluctuation while identifying the point of convergence. The procedure is as follows: 1) Eliminate the data point xj,{xj, xj≥xmax-d×j} from the precipitation time series {xi, i = 1, 2,..., n}, where xmax is the maximum, xave is the average, and d is the interval space (we set d =1 mm in this case). 160 2) Compute the selected parameters.
Using the obtained parameter series, we applied the segmentation algorithm proposed by Bernaola Galvá n et al. (2001) to determine the point of abrupt change, which we define as the EPT. The segmentation algorithm is based on the calculation of the statistic t of each data point in a series or subseries: 165 where sD = [(sL 2 + sR 2 ) / (NL + NR -2)] 1/2 (1 / NL + 1 / NR) 1/2 is the pooled variance. The μL/μR, sL/sR, and NL/NR are the mean, standard deviation, and number of points from the data to the left/right of the series, respectively. The significance level P() of the largest value tmax obtained from Eq. (11) is defined as the probability of obtaining a value equal to or less than τ within a random sequence (Swendsen and Wang, 1987) If the significance exceeds a selected threshold P0 (usually taken to be 95%), an abrupt point is selected and the series is cut at this point into two subsets.
The pooled variances and the abrupt points of the α and γs series are shown in Figures 1b-1d. The abrupt point, where sD differs from its left-side points but is approximately equal to its right-side points, is selected to be the EPT. As shown in Figure 1, for 175 the determined EPT, there are a deep slope of the pooled variance on the left while a gentle slope on the right of daily precipitation 37.1 mm. Thus, the EPT for Xingxian station is estimated as 37.1 mm/d, and 90 EP events have occurred over the past 55 years ( Figure 1c).

EP indices
All of the variables characterizing spatiotemporal EP over the LP are shown in Table 1. For individual station, annual 180 precipitation at each station was accumulated using daily data. Annual EP is accumulated using daily EP determined by EPT, and EP frequency (EPF) is the number of daily EP events. Mean annual EP (MEP) is averaged from annual EP interpolated using ArcGis 10.1. For a year with n EP events, EP intensity (EPI) is calculated by the equation where EPi represents the magnitude of EP event i, respectively. 185 As noted by IPCC (2007), neither the cases of high frequency with low intensity, or low frequency with high intensity, can reflect the severity of EP events given a long-term time series for an area, while the severity of EP (EPS) events relies on both intensity and frequency. Consequently, we examine the extreme precipitation intensity (EPI), extreme precipitation frequency 8 (EPF, defined below) and EPS to characterize the spatiotemporal nature of EP over the LP. To obtain the concordant EPS, each annual EPF/EPI series is standardized to the range 0 to 1 using equation (10): 190 where Xmin and Xmax represent the lowest and highest annual EP frequency/intensity, respectively. The annual EPS for each station (0 ≤ EPS ≤ 1) is calculated from the standardized EPI and EPF by where k1 and k2 are the weight coefficients of frequency and intensity influencing EP severity, respectively, and k1 + k2 = 1. In 195 this paper, k1 and k2 are set 0.5 according to IPCC (2007).
In addition, we compute the long-term mean daily precipitation (MDP) and the long-term accumulated daily EP events (ADEP) to help characterize the intra-annual precipitation and EP. As shown in Table 1, MDP is averaged from all 87 stations, and ADEP is accumulated from 87 stations over the past 55 years.

Spatiotemporal EP presentation 200
The spatial distributions of the EP indices were derived by interpolation via Kriging (Oliver and Webster, 1990), using data observed at the gauging stations. All spatial analyses were carried out using the ArcGis 10.1 software. Spatiotemporal trends for annual EP indices were computed for each pixel using the least squares method. Following Stow et al. (2003), the trend is where m is the total of years, Pi is value of the pixel in the j-th year.

Study area
The LP (640 000 km 2 ) is a typical arid and semiarid area located in the middle reaches of the Yellow River (750 000 km 2 ) and 210 characterized by a continental monsoon climate. Elevations range from 84 m to 5207 m ( Figure 2). The desert-steppe, typical steppe, and forest steppe (deciduous broad-leaf forest) zones are distributed from northwest to southeast, and correspond to mean annual isohyets of 250, 450, and 550 mm in the arid, semiarid, and semi-humid areas, respectively. The continuous loess covering ranges from 100 m to 300 m in thickness on the mountains, hills, basins, and alluvial plains of different heights. The northwestern part of the region is dominated by flat sandy areas. The middle and southeastern parts are characterized by EP 215 induced water-erosion landform (Zhang et al., 1997), with a rugged undulating ground surface that is broken, barren, and dissected by gullies and ravines (Cai, 2001). EP-induced flooding episodes occur occasionally in the summer, with sediment concentrations generally exceeding 300 kg/m 3 and have been observed to be as large as 1,240 kg/m 3 . The hyper-concentrated flooding has historically resulted in numerous disasters, with severe consequences to people and livestock (Zhang et al., 2017).
The amount of soil erosion has been estimated to be larger than 2000 to 3000 million tons per year (Tang, 1990). Soil erosion 220 has resulted in the density of gullies and ravines in the LP being larger than 3-4 km/km 2 , with the maximum exceeding 10 km/km 2 .

Database
To conduct the EP assessment, we used daily data available for 87 national meteorological stations in and around the LP  (Kalnay et al., 1996) were also used in this study. The variables selected for analysis were the monthly mean geopotential height, monthly mean wind, daily mean sea level pressure and daily 235 mean wind from 1961 to 2013 on a 2.5°×2.5° spatial grid (http://www.esrl.noaa.gov/psd/). 3b); these EPT isohyets are generally consistent with these of mean annual precipitation. Figure 3b indicates that the area around the Dongsheng station is a regional EP center since that the EPTs around the station are higher than these of the surrounding areas, whereas the isohyets of mean annual precipitation are smooth. The spatial distribution of MEPs is also similar to that of mean annual precipitation, increasing from northwest to southeast and ranging from 35 to 138 mm· yr -1 , as shown in Figure 3c. Maximums of MEP occur in the southern and southeastern LP. 245 Figure 3d indicates that EPFs over the region during the past five decades have ranged from 54 to 115 d (i.e., mean annual EPF ranged from 1.0 to 2.1 d). Notable occurrences of high EPF can be seen in and around the Ziwuling Mountains in the mid-southern LP, while the highest frequency occurred at the east of the Fenhe Valley in the southeastern LP. Meanwhile, the lowest frequency occurred in the northwestern LP, the western Muus sandy land, and the western Liupan Mountains. This high value of EPI in part explains why this area is characterized by very serious soil erosion that releases more than 2billion ton sediment annually into channels of the Yellow River. Figure 3f indicates that the spatial distribution of EPS in the LP increased from northwest to southeast, ranging from 0.27 to 0.66, but with the highest EPSs centered in the southeast Central LP. The areas with highest EPSs covered the basins of the 260 Jing, Luo, and Fen rivers. Although an EP event always occurred over a small range, the spatial maps of EPI, EPF, and EPS indicate that the areas with serious EP events are regularly distributed.

Spatiotemporal variation of EP
Our results (Figure 4a) indicate that 91.4% of the LP was characterized by a negative annual precipitation trend over the study period, whereas only 8.6% of the total area presented a positive trend. There were 9 out of 87 stations showed significant 265 negative trend while 2 stations showed positive trends (p < 0.1). At the same time, the spatiotemporal trends of the annual EP ranged from -0.78 to +0.48 mm· yr -1 (Figure 4b), with 23.8% of the total area showing a positive trend, with increased annual EP distributed mainly in the southwestern LP (west of Lanzhou) and the mid-southern LP (Beiluo and Jing river basins and an area around the Xingxian station). Meanwhile, the annual EPF changed by -0.6 to +0.5 d over the past 55 years, with a change rate ranging from -1.2 × 10 -2 to +0.95 × 10 -2 d· yr -1 , as shown in Figure 4c. Of the 87 stations, 4 stations showed significant 270 negative trend of EPF while 3 stations showed positive trends (p < 0.1). The areas with a negatively trending EPF covered 86.4% of the total area while the areas with positively trending EPF covered 13.6%, the latter occurring mainly in the southwestern LP (around the Xining station) and in the areas around the Xi'an and Xingxian stations. The areas with notably decreasing trends occurred mainly in the mid-west and southeast regions of the LP. Figure 4d indicates that the changes of annual EPI ranged from -0.18 to +0.27, with a changing rate ranging from -0.34 × 275 10 -2 to 0.52 × 10 -2 yr -1 . We found that 34 of the 87 stations showed upward slopes (5 stations with a significance level p < 0.1), and 53 stations (4 stations with a significance level p < 0.1) showed negative slopes. As shown in Figure 4d, areas with positive trends of EPI accounted for 42.2% of the total area, with the areas delineating by the Wulate -Yan'an-Huashan stations and the Jingtai-Xiji-Tianshui stations, as well as the area west of the Minhe station. The areas with a negative slope covered 57.8% of the total area. 280 Figure 4e indicates that the annual EPSs changed by -0.09 to +0.07 during the study period, with rates varying from -0.34 × 10 -2 to 0.52 × 10 -2 yr -1 . Of the 87 stations, 39 stations showed a positive slope (3 stations with a significance level of p < 0.05), while 54 stations exhibited a negative slope (2 station with a significance level of p < 0.05). The areas with increased EPSs covered 25.4% of the total area and were mainly found in an area delineated by the Wuqi, Tianshui, and Huashan stations and an area west of the Xiji station. The areas with negative trends accounted for 74.6% of the total area. 285 The trend estimates computed for annual EP, EPF, and EPI are associated with strong uncertainty. For instance, the upward trend of annual EP in and around the Xingxian station relied heavily on the upward trend of the EPF and not the downward trend of the EPI. The EPF around the Changwu station decreased, but both the annual EP and EPS increased with the upward trend of the EPI. However, nearly all the areas with positive trends for annual EP, EPI, EPF and EPS had a negative trend in annual precipitation (Figure 4). It should be noted that 62.1% of the LP area with negative trend in annual precipitation has 290 more than one EP indices with positive trends, potentially indicating the risk of more serious hazardous situations. Precipitation from June to September accounts for 72% of the total amount, while 91% of the total EP events occur from June 295 to August. According to the fitted curve ( Figure 5), the highest MDP occurred on July 26, which is 11 days earlier than the maximum ADEP on 6 August. Based on fitting the four-parameter Weibull curve (p < 0.0001), the MDP for the 224 days from March 26 to November 4 accounted for 95% of the mean annual precipitation. Meanwhile, the ADEP from May 21 to September 18 accounted for 95% of the total EPF. Therefore, high concentration of amount of daily precipitation into a limited period results in a significant alternation of wet 300 and dry seasons in the LP. In addition, low precipitation but with annual alteration of dry and wet seasons, and highly concentrated intra-annual EP events with an occurrence 11 days earlier that the wettest days, contributes to a fragile ecoenvironment subject to severe natural hazards. Specifically, lower precipitation but highest EPI and EPS are responsible for the most severe hazard situations in the Central LP, such as soil erosion.

Atmospheric circulation factors for the spatial variation of extreme precipitation 305
Atmospheric circulation is the leading factor causing the above phenomena. The LP is located in the East Asian monsoon region. According to the average sea level pressure and winds in winter from 1961 to 2015 (see Figure 6a), the dry winter in the region is influenced by the interactions between two high pressure areas in Southwest China (the Tibet Plateau high pressure system) and North China (the Mongolia high pressure system). The prevailing East Asian winter monsoon (which has a northnorthwest direction) circulates in East China and brings cold and dry airstreams. In contrast, the summer climate of the LP is 310 affected by interactions between two high pressure systems, the Pacific high pressure and Tibet Plateau high pressure systems. Figure 6b shows that the prevailing East Asian summer monsoon (which has a south-southeast wind direction) brings warm and humid maritime airstreams that spread from the West Pacific to Central China. However, the Tibet Plateau high pressure has a notable effect on the climate of the northwestern LP, and the airstream humidity decreases gradually as the distance from the Pacific increases. The resulting effect and decreased humidity is to form a vast arid region in Northwest China, including 315 the northwestern LP, with a prevailing wind direction of west-southwest. This explains why precipitation decreases from southeast to northwest, and precipitation is scarce in the northwest LP.
Nevertheless, tropical cyclones occasionally enter the central LP, accompanied by EP events. For instance, in August 1996 a Western Pacific cyclone landed in the southeastern coastal area of China and weakened gradually as it moved northwest, as shown by the 1000-hPa geopotential height and winds in Figure 6c. Plenty of rainstorms or intense rainfall events accompanied 320 the cyclones occurred in its transit area. On 3 August 1996, the weakened cyclone reached the southeastern LP, as shown in Figure 6d. However, under the control of the Tibet Plateau high pressure, the central LP is generally the northwestern boundary to which the tropical cyclone can reach. As shown in Figure 6d, the cyclone was blocked from entering the northwestern LP, moved towards the northeast, and gradually dissipated. These phenomena illustrate why this region has limited precipitation but severe EP events. 325

Rationality of Spatial EP Characteristics
Natural hazards related to EP can be divided into two categories: (1) hazards accompanied by EP, and (2) hazards that follow the occurrence of EP. For the former, one focus is the dependence of EP and storm surges in the coastal zone. Using such dependence structures, EP and storm surge can be quantified to provide information for successful hazard management 330 (Svensson and Jones, 2004;Zheng et al., 2013). For the latter, the LP is such that the area suffers from EP-induced natural hazards that exceed the general tolerance of the natural environment, existing ecosystems, human life and social economy. In this case, the rational characteristics of EP responsible for spatial hazards can be studied.
Here, we use the widely distributed soil erosion to verify the rationality of our results. According to the universal soil loss equation (Wischmeier, 1976), the rational characteristics of EP should correlate well with soil erosion and vegetation coverage 335 13 (Figure 7). To examine this, partial correlation analyses were performed between soil erosion intensity and EPI/EPS, and with the vegetation coverage. Our results indicate that water-based erosion intensity correlates significantly with vegetation coverage (negatively, Figure 7a) and EPS (positively, Figure 3f); the related coefficients are -0.61 (p < 0.001) and 0.53 (p < 0.001), respectively. For the correlations between water erosion intensity and vegetation coverage, and with EPI (Figure 3e), the coefficients are -0.58 (p < 0.001) and 0.76 (p < 0.001), respectively. This finding demonstrates the rationality of our results. 340 Note that, the higher correlation between EPI and soil erosion agrees with the results of plot experiments by Tang (1993), who noted that high-intensity precipitation is the primary driving force of erosion. Zhou and Wang (1992) divided the LP into three zones of raindrop kinetic energy (<1000, 1000-1500, and 1500-2000 J· m -2 · yr -1 , respectively), based on observations of the raindrop kinetic energies of rainstorms during 1980s. We found that the 30 and 35 mm/d EPT contours closely overlap with the raindrop kinetic energy contours of 1000 and 1500 J· m -2 · yr -1 . Further, 345 soil erosion in the LP in recent decades has been found to be approximately 5000-10 000 t· km -2 · yr -1 (Ludwig and Probst, 1998;Shi and Shao, 2000). Such high rates of sediment erosion are generally induced by several rainstorm events during the year, with the top 5 daily sediment yields accounting for 70%-90% of the annual total soil loss (Rustomji et al., 2008;Zhang et al., 2017). For instance, a 200-year precipitation event in Wuqi on 30 August 1994 induced a flooding event with a daily sediment concentration of 1060 kg/m 3 . The streamflow was 2.41 times of the mean annual streamflow from 2002 to 2011, and the 350 sediment load was equivalent to 9.6% of the total sediment yields from 1963 to 2011 (Zhang et al., 2016). Therefore, the EPF obtained in this study, about twice every year on average, is rational to explain such serious sediment erosion. In the LP, spring drought is the limit factor for vegetation (especially herbaceous vegetation) recovery from winter every year, and grass generally germinates on an extensive scale after the first effective rainfall event (>5 mm) in spring (Cai, 2001;Tang, 1993;Tang, 2004). However, as shown in Figure 5, the highest EPF occurred 11 days earlier than the day of maximum daily 355 precipitation in the LP. This means that, the days on which the LP experiences most serious EP events, tend to be days when precipitation is less. In other words, every year, the vegetation has not sufficiently recovered when frequently EP events occur in the LP. Such an intra-annual distribution of precipitation is one of the climatic reasons why there is serious soil erosion in the semiarid LP. Further, the sparse spatial nature of precipitation is insufficient for the growth of high-coverage vegetation, especially in the northwestern LP ( Figure 7a). However, the highest EPI provides the strongest erosion force, which contributes 360 to the severe rates of erosion (Figure 7b) in the Central LP. These results of EP responsible for hazardous situations in both spatial and temporal are important for sustainable catchment management, ecosystem restoration, and water resources planning and management within the LP. Given that 62.1% of the total LP with negative trend of annual precipitation has one or more positive EP indices, the underlying upward trends of water erosion and sediment yield should be taken into account in catchment management efforts. 365 14

Uncertainty in EP identification
The uncertainties in identification and assessment of EP events come from two aspects: (1) the stochasticity in climate (Miao et al., 2018) and (2) the methodology (Papalexiou et al., 2013). For the former, significant spatiotemporal variations occur in EP events as a result of varying geographical and meteorological conditions (Pinya et al., 2015). Extreme precipitation events are relatively rare, poorly predictable, and often with short duration, thus resulting in uncertainty in EP event identification. In 370 the method section, the uncertainties in EPTs determination from parametric and non-parametric methods were discussed. Parametric methods require a predetermined threshold value, above which the data can be chosen as the EP series if the data series passed the goodness-of-fit test. As shown in Figure 9, both fixed values and percentiles were adopted to preset EPT. The selected rainfall series data were fitted to the gamma, GPA (generalized Pareto distribution), Gumbel, Pearson type III, GEV 385 (generalized extreme value distribution), and the GNO (generalized normal) distributions, whose parameters were estimated with the L-moments method (Haddad et al., 2011) at a 0.05 significance level, using goodness-of-fit tests including K-S (Kolmogorov-Smirnov test), A-D (Anderson-Darling K-Sample test) and C-S (Pearson's Chi-squared test) tests. As shown in Figures 9a1-9a3 and 9b1-9b3, the results of the three tests are similar but different in details for the preset fixed value and percentile thresholds. 390 Further, these results for different distribution functions are quite different from each other. As shown in Figures 9a1-9a2, by the K-S and A-D tests, the passing rates from GEV, Gumbell, GNO and pearson type III distribution functions are high while there is very low passing rate from GPA and Gamma functions. In addition, the passing rates are different between preset methods of percentiles and fixed values. As shown in Figure 9, the GEV and Gumbell distribution function have high passing rates for EP series obtained by preset percentiles when percentile ≤ 99% (Figures 9a1-9a2), whereas the passing rate for these 395 series obtained by fixed values decrease with increasing values (Figures 9b1-9b2). We also found that these distribution functions are not sensitive to the percentile or fixed value changes. These findings indicates that the fitting accuracy can be greatly affected by the selection of the extreme value distribution functions, goodness-of-fit tests and methods for EPT preset. Liu et al. (2013) noted that the fitting accuracy is also affected by the size of rainfall series. So, unavoidably, applications of parametric methods also depend on personal subjectivity and empiricism. We have tried to explore EP using fixed values in 400 spatiotemporal variation of precipitation analysis in the LP; however, we found that the results could not well explain the rainfall-induced natural hazards or agree with these plot experiments in spatial (Wan et al., 2014), and the same to these results obtained by fixed values (Xin et al., 2009) and percentiles (Li et al., 2010a;2010b). These uncertainties may the reason why prior studies of EP over the LP tend to disagree with each other.
As noted by Pandey et al. (1998) and Douglas and Barros (2003), these methodological uncertainties arise due to the wide 405 gap between mathematical modelling and the physical understanding of precipitation processes. As previously mentioned, the multifractal technique can be used to describe the statistical probability and physical processes associated with observed data (Lovejoy and Schertzer, 2013;Tessier et al., 1996), while the scale invariance of multifractals enables the multifractal technique to also overcome the influence of the sample size (Pandey et al., 1998;Tessier et al., 1996). Further, the segmentation algorithm helps to overcome the problem of uncertainty. In the present study, the general correspondence and the specific 410 divergences between EPT and precipitation isohyets (Figure 3) further exhibits the varying meteorological and geographical influences. As shown in Table 2, by fitting EP series derived by universal multifractals to the six distribution functions, the 100% passing rate of goodness-fit test strongly supports that universal multifractal approach is advanced in identifying EP events. Overall, the universal multifractal method provides a much superior approach to addressing uncertainties and providing a unique set of EPTs. 415

Conclusions
We have proposed an approach that integrates universal multifractals with a segmentation algorithm to enable identification of EP events, and thereby to assess its spatiotemporal EP characteristics in the LP, using data from 87 meteorological stations from 1961 to 2015. We find that the spatial distribution of the EPTs increased from 17.3 mm/d in the northwestern to 50.3 mm/d in the southeastern LP. Similarly, the MEP increased from 35 mm to 138 mm/yr, with the maximum MEP occurred in 420 the southern and southeastern LP. The EPF over the LP was within a range of 54-116 days over the last 55 years. Notable occurrences of EPFs mainly observed in the mid-southern and southeastern LP. An examination of atmosphere circulation patterns demonstrates that the central LP is the boundary where tropical cyclones enter the inland China, resulting in the highest EP intensity and EP severity in this area. Correlation analysis significantly supported the reasonability of the spatial estimates of EP characteristics that are responsible for hazardous situations over the LP. The climate factors for the most serious 425 hazardous situations in the LP especially in the Central LP come from the low precipitation, the highest EPI and the highly ADEP concentrated 11 days earlier than the wet season.
Spatiotemporally, annual EP increased in the southwestern and mid-southern LP. The areas with a positive EPF trend occurred in the southwestern LP and the areas around the Xi'an and Xingxian stations, whilst the areas with a positive trend of EPI among the Wulate-Yulin-Yan'an-Huashan stations and the Jingtai-Xiji-Tianshui stations, as well as the area west of 430 the Minhe station. The annual EPSs with increased slope covered an area delineated by the Wuqi, Tianshui, and Huashan stations and an area west of the Xiji station. Overall, the areas with upward trends of the annual EP, EPF, EPI, and EPS accounted for 23.8%, 13.6%, 42.2%, and 25.4% of the LP area, respectively. It should be noted that 62.1% of the LP area with negative annual precipitation experienced upward trends of one or more EP variables. It can be concluded that EP over the LP intensified, potentially imposing a risk of more serious hazardous situation. Sustainable countermeasures should be considered 435 in the catchment management to address the underlying hazards.
In conclusion, the universal multifractal approach considers both the physical processes and their probability distribution, and thereby provides an approach to overcome uncertainties and identify EP events without the need for empirical adjustments.
This approach is thus useful for application to spatiotemporal EP assessment at regional scale.