Contribution Analysis of the Streamﬂow Changes in Selected Catchments on the Loess Plateau, China, Using Multiple Budyko-Based Approaches

: A better understanding of how streamﬂow interacts with climate change and human activities would contribute to the efﬁciency and effectiveness of water resources management. Speciﬁcally, quantifying the climate and human contributions has widely been used when attributing streamﬂow changes. However, only a few previous studies compared the results derived by different methods that are currently available, and even fewer studies have ever had a close look at the uncertainties induced by various estimations of evapotranspiration. This research ﬁrst examined the streamﬂow changes for 12 catchments on the Loess Plateau in China during the period of 1961–2018 with Mann–Kendall test and relevant statistical measurements. Then, 8 Budyko-based climate elasticity methods, each with 13 estimations of evapotranspiration, were used to quantifying human and climate contributions to streamﬂow change in the study area (i.e., 104 pairs of values for human and climate contributions for one catchment). The results showed that signiﬁcant declining trends could be found in 11 of the 12 catchments studied. In terms of contribution rates, human activity has been shown as the major contributor to the streamﬂow decrease (60–90%) compared to climate change (10–50%). By comparing the contribution results derived by possible combinations of attribution method and evapotranspiration estimation, the variability due to different Budyko-based methods being used seems to be related to geographical location and climate. Although the spatial pattern of variability due to different estimations of evapotranspiration is not obvious, it is necessary to consider the uncertainties induced when launching contribution analysis over speciﬁc regions.


Introduction
Water resources are closely related to agricultural production and human well-being. Meanwhile, both climatic changes and human interferences can influence the availability of water resources [1]. Climate change has been observed to have impacts on the global water cycle and alter the patterns of regional extremes (e.g., flood and drought) [2]. On the other hand, human activities also affect water dynamics through the water extraction and production activities that lead to changes in land use and land cover [3,4]. Climatic and anthropogenic factors sometimes come into play together, which leads to intensified impacts on natural flows and extreme conditions, such as drought [5]. In this manner, it is necessary to better understand the impacts of either internal or external factors on water resources and the underlying mechanisms to help decision makers manage water resources in a more sustainable way. Streamflow, also known as runoff, is one of the key components in hydrologic cycles. It is important to investigate the influence of environmental changes (due to either climatic or anthropogenic factors) on streamflow to assess how regional water resources would respond to the changes and to generate adaptation strategies for the changes [6]. Previous studies have investigated the changes in volume of runoff and its flow regime, aiming at identifying the characteristics of streamflow changes and the factors driving streamflow changes [7][8][9]. There are numerous approaches to describe the variations and characteristics of streamflow, such as indicators of hydrologic alteration (IHA), range of variability approach (RVA), and flow duration curve (FDC)-or hydrographbased hydrologic alteration assessments [10][11][12]. Over a long period, both climate change and human activities contribute to streamflow changes [13]. Therefore, it is vital to separate impacts by different factors when analyzing the causal relationship between streamflow changes and influencing factors.
Previous studies have applied several different methods to partition climate and human contributions. Specifically, the partitioning methods that are currently available can be categorized into physically based hydrological models, statistical methods, and climate elasticity methods based on the Budyko hypothesis [13,14]. Li et all (2009) claimed a predominant impact of land use changes on hydrology within a catchment located at LP by applying SWAT model in their research [15]. Alike models are often affected by the lack of high-resolution datasets and the uncertainties raised by parameterization and model structures [16,17]. By using linear regression methods, Xu (2011) reported a higher contribution of climate change to the variability of streamflow [18], while other studies showed a relatively lower contribution of climate change, in which different study periods and research methods were applied [15,19]. Recently, Budyko-based climate elasticity methods have been increasingly used in the assessment of streamflow changes, because they are physically centered and only involve simple calculations [20,21].
In terms of the Budyko framework, previous studies have proposed several Budykobased hypotheses that define the annual water balance (i.e., the ratio of actual evapotranspiration to precipitation) as a function of the ratio of potential evapotranspiration to precipitation (i.e., dryness index, ϕ) [22]. The climate elasticity methods based on those hypotheses can then be used to quantitatively separate the impacts of climate change and human activities on hydrology. However, few previous studies have investigated streamflow changes by quantifying human and climate contributions with multiple Budyko-based methods and assessing the uncertainties involved in the results of quantitative contribution [13]. One of the input parameters for the Budyko-based climate elasticity methods is potential evapotranspiration (E 0 ) the value of which, unlike other parameters, can be obtained from gauging stations; the E 0 being used here is often estimated using radiationbased methods, such as the Penman-Monteith (P-M for short) method recommended by Food and Agriculture Organization of the United Nations (FAO). Different estimations of E 0 can lead to uncertainties in results but there is a lack of studies focusing on the comparison contribution analyses using different Budyko-based methods and varying inputs of E 0 .
The Loess Plateau (LP) in China, due to its long history of agricultural practice, low level of vegetative coverage, and the regional climate characterized by intense rainfall in summer, has high soil erosion rates and heavy sediment loads, which, in turn, make the ecosystem here extremely vulnerable to climate change [23]. Numerous studies have investigated the historical climate change using observational datasets [24,25]. For LP, Sun et all (2015) examined the mean and extreme values of temperature and precipitation from 1961 to 2011 and found that there was a general warming trend across the plateau while the magnitude varied spatially; on the other hand, more than one third of the plateau had witnessed a decrease in precipitation, which can be attributed to decreases in the frequency of the event or the intensity of the event [26]. Given the background of global warming, the temperature of LP increased at an annual rate of 0.02 • C since the 1950s, and along with the temporal inequality in precipitation, it contributed to the drying trends identified in almost half the plateau [27,28]. From the 1950s to the 1980s, efforts have been taken to increase the resilience of the LP ecosystem, including a series of water-soil conservation practices (e.g., afforestation and grassing to alter land surface conditions and micro-topography) and hydrological engineering events, such as dam constructions [29,30]. Then, in the 1990s, after experiencing severe drought and subsequent floods, two major conservation actions were launched, aiming at controlling soil erosion and sediment transport, known as the Natural Forest Conservation Program (NFCP) and the Grain-for-Green Program (GFGP) [31]. Those nationwide actions promoted the existing conservation measures, such as building check-dams, converting farmland to forest or grassland, and making investments on ecosystem services, termed as the Ecological restoration (ER) measures. Those ER activities have been proven to be effective in vegetation recovery and reducing sediment [32,33]. The rainfall interception and actual evapotranspiration in this region increased consequently. Meanwhile, an obvious decreasing trend in the streamflow within LP has been noticed, which affects the ecosystem in a negative way [34][35][36]. Thus, relevant studies are typically interested in the water resource changes within LP and the underlying mechanisms [21,37,38].
Considering the ongoing attempts, or those to be raised in the future, to protect the fragile ecosystem of LP, a better understanding of the impacts of human interference and climate change on the streamflow change would be essential. Therefore, this research employed eight Budyko-based climate elasticity methods, each with 13 estimations of E 0 , to analyze the contributions of climate change and human activities to the streamflow changes in 12 catchments on LP during the period 1961-2018. The objectives are to (1) have more insights into the variabilities of LP streamflow; (2) capture the characteristics of the interactions between climate, human activity, and streamflow change within the study area; (3) discuss the uncertainties in the contribution analysis results that associated with different methods and parameter values.

Study Area and Data
The Loess Plateau (approximately 35-41 • N, 102-114 • E) is located in the middle reaches of the Yellow River basin in North China. Approximately 8.5% of the Chinese population live on this plateau, which covers an area of 632,520 km 2 (approximately 6.6% of the entire land area of China). The plateau is dominated by arid and semiarid climates with relatively higher evaporation rates, while it is affected by both East and South Asian monsoons that mostly bring rainfall during the period from June to September. The highintensity rainstorms, along with the coverage of highly erodible loess, accelerate the soil erosion process, especially for the regions with slopes that exceed 10 • [21,39,40]. This study will investigate the hydrology in 12 catchments ( Figure 1) across the regions with varying precipitation, land cover types, and land surface gradients [13,23].
For the study area, monthly streamflow data from 1961 to 2018 were acquired from the Yellow River Conservancy Commission (YRCC). Daily precipitation data was provided by the China Meteorological Data Service Centre (CMDC). The potential evapotranspiration was calculated by P-M method and 6 other radiation-base methods, each of which was processed with two different parameterizations (more details are discussed in Section 3.2.) [41,42]. ArcGIS was applied to conduct linear interpolation to spatially average meteorological variables across the whole region of interest. In addition, the interannual variabilities in normalized difference vegetation index (NDVI) with respect to individual catchments were checked to examine the correlation between vegetation cover and streamflow changes. The data were obtained from MODIS with a spatial resolution of 1 km.

Mann-Kendall Test
The Mann-Kendall (M-K) test has been widely used to detect trends in the time series of hydrological variables [8,43,44] as well as vegetation indices [45,46]. This test is based on a test statistic S, defined as: where x i and x j are observations, n is the sample size, and sgn is the sign function. When n is larger than 8, S can be approximated as the normal distribution with a mean value of zero [47]. The variance of S can be calculated as: where t i is the number of ties of extent i. Then, the standard test statistic Z can be calculated: By conducting the M-K test, there will be an increasing trend when Z is greater than zero and a decreasing trend when Z is less than zero.

Indicators of Hydrological Alteration (IHA) and Range of Variability Approach (RVA)
IHA has been widely used to measure streamflow changes and, based on that RVA, was developed to have more insights into this hydrological variability. IHA includes several groups of statistics describing different properties of stream flow changes, such as the magnitude of monthly water conditions and the magnitude and duration of annual extreme water conditions [48]. Here, the first group of statistics was selected considering the temporal resolution of data used in this research. Specific parameters included in the selected group are listed in Table 1. Table 1. IHA statistics (Group 1). For the baseline period of a long-term time series of hydrological data, RVA would typically set the range of IHA variations as the interval between the 25th and 75th percentile, or the interval defined by one standard deviation from the mean IHA [10]. Then, the IHAs for the evaluation period will be examined to see if they fall into the pre-defined range. Those IHAs falling outside the range reflect the alteration, which can be measure by the alteration ratio as:

IHA Statistics Hydrological Parameters
where N o,i is the number of years when the stream flow is within the range, the subscript i denotes the corresponding IHA statistics, and the Ne is the expected number of years when the stream flow is within the range, the value of which is typically half the total number of years. Furthermore, the alteration ratios for the high season, low season, and total year can be computed as: The absolute value of D between 0-33% is deemed as mild alteration, between 33-67% as moderate alteration, and between 67-100% as severe alteration.

Budyko Hypothesis
On a long-term average, the Budyko hypothesis ignores the change in water storage of a catchment and uses physical principles to show that the water balance (represented by the ratio E/P, where E is actual evapotranspiration and P is precipitation) is largely controlled by the water-energy balance (represented by the ratio of E 0 /P). Various equations, known as Budyko-type equations, have been proposed along the years as the solution to the Budyko hypothesis. Wu et all (2017) tabulated a total of 8 mathematical functions that represent the Budyko hypothesis, respectively (Table 2) [13]. Some of these equations are non-parametric, while others include parameters that are connected with catchment characteristics, such as soil type, topography, and vegetation.  Climate elasticity method was first proposed to identify the relationship between streamflow (Q) changes and their influencing factors (ε X = ∂Q/Q ∂X/X ) [49]. Based on the Budyko hypothesis, numerous studies have attempted to explain the variability of streamflow by assuming that the changes in E 0 , P, and human activities are independent throughout the catchment of interest [50,51]. In this manner, the streamflow changes induced by climate change (Q c ) can be expressed as:

Function Identifier f(∅) Parameters
According to the concept of climate elasticity, Equation (9) can be rewritten as: Then, the P elasticity of Q, 0, and the E 0 elasticity of Q, ε E0 , can be given as: Following this, the climate-induced changes in streamflow can be calculated as: It is worth noting that all the precipitation, potential evapotranspiration, and streamflow in the equation above are annually averaged over a long period.

Decomposition Method
As suggested by Wu et all (2017) and Liang et all (2015) the study period  can be split into two sub-periods: the baseline period that was less disturbed by human activities (1961)(1962)(1963)(1964)(1965)(1966)(1967)(1968)(1969)(1970) and the evaluation period (1971-2018), within which human activities were greatly intensified [13,23]. In this research, the difference between the mean annual streamflow in the baseline period and that in the evaluation period is regarded as the total change in the streamflow (∆Q). Thus, the relationship between the two sources that contribute to the streamflow changes can be expressed as: where ∆Q h is the streamflow changes caused by human activities. Furthermore, the contributions of climate change (C c ) and human activities (C h ) can be given, in percentage, as:

E 0 Estimation Methods
Potential evapotranspiration (E 0 ) is an important process in both water balance and energy balance. The assessment of a climate's influence on water resources could be associated with E 0 estimation. However, E 0 is controlled by many factors, including humidity, wind, and radiation fluxes, which may lead to uncertainties in the results of estimation models. While there are several estimation methods available, most of them were developed focusing on different regions or climate types and were based on specific assumptions [52][53][54][55][56][57][58][59]. Therefore, it is necessary to compare the estimation efficiency in various circumstances to determine the feasibility of individual methods. In this study, we would like to examine the sensitivity of contribution analysis to the choice of E 0 estimation method. More details about the methods involved in this study can be found in Table 3. Table 3. E 0 estimation methods.

Method
Function

Changes in Hydrometeorological Variables
The M-K test was applied to examine the trends in the time series of three hydrometeorological variables (i.e., Q, P, E 0 ) for the 12 catchments at both interannual and intraannual scales. Specifically, the intra-annual analysis includes trends for high-flow (June to September) and low-flow seasons (October to May), respectively. As shown in Table 4, no significant trend can be captured for P over the study area; in all the catchments, E 0 increased significantly in the low-flow seasons at rates ranging from 2.1 mm/yr (Kuye) to 2.5 mm/yr (Beiluo), while only the E 0 in Xinshui shows a significant upward trend in high-flow season. In terms of annual E 0 , again, not all catchments exhibit significant upward trends. On the other hand, streamflow shows downward trends in both high-flow and low-flow seasons in most catchments, except that no such trend was captured for Dali and Qingjian. Furthermore, Dali turned out to be the only catchment without a significant downward trend found at the interannual scale.
In this table, ↓ denotes that the Z-statistic of Mann-Kendall test is smaller than zero and indicates a downward trend; ↑ denotes that the Z-statistic of Mann-Kendall test is larger than zero and indicates an upward trend;indicates that there is no significant trend detected; *** indicates that p value is less than 0.001, ** less than 0.01, and * less than 0.05; for the RVA results of streamflow, L denotes the mild alteration, M denotes the moderate alteration, and H denotes the severe alteration.
RVA was then applied to the streamflow time series of all 12 catchments for further investigation. The results are also shown in Table 4. In general, the annual streamflow in the catchment Fenhe and Jialu were severely altered; those in the catchment Jinghe and Beiluo were mildly altered; the rest of the catchments have witnessed moderate alternations in the streamflow. At a finer temporal resolution, alterations are more noticeable in high-flow seasons than in low-flow seasons. Specifically, in high-flow seasons, three catchments in the study area have their streamflow severely altered, and the streamflow of the other catchments were also moderately altered. In comparison, more than half the catchments have their streamflow only mildly altered in low-flow seasons. The temporal patterns of streamflow changes in the catchments Fenhe and Jialu are so different from the others' that most server alterations to the streamflow in those two catchments appeared in the low-flow seasons. The time series of streamflow more explicitly depicted the changes (Appendix A, Figure A1). In general, the streamflow within most catchments studied in this research have fluctuated along the years, but the overall trends are downwards. Compared to low-flow seasons, the streamflow fluctuated more obviously in high-flow seasons, especially for the catchment Gushan.

Contribution Analysis of Streamflow Changes
In this section, the quantitative contributions of climate change and human activity to streamflow changes in LP are analyzed by comparing results derived by different Budykobased methods as well as various estimations of evapotranspiration (E 0 ). This research applied eight different Budyko-based methods and 13 different estimation methods for E 0 (except for the FAO method, all other methods have corresponding versions that have been tuned by adjusting specific parameters involved). In other words, there are a total of 104 pairs of values for the contribution of climate change and human activity to stream change in each of the 12 catchments being investigated in this research (quantitative attribution of stream flow change for each catchment is not shown here). Figure 2a  To control the variables involved in analysis, the results with different E 0 estimation methods involved will be averaged (i.e., there is only one pair of human and climate contributions left) when comparing various Budyko-based methods, and similarly, the results derived by different Budyko-based methods will also be averaged when comparing various E 0 estimation methods.

Quantitative Attribution Using Various Budyko-Based Methods
As shown by Figure 3, for a specific Budyko-based method, the results were averaged over calculations using different E 0 . Note that the sum of human and climate contributions equals 100% by design. A red-white-blue color ramp was applied to the results of human contribution, and the same color ramp was reversed for the results of climate contribution. Consequently, in both upper and lower panel, the color that is closer to the red side of the color ramp indicates that human activity contributes more to streamflow change compared to climate change. From Figure 3, we can tell that, for catchments such as Tuwei, Gushan, and Kuye, the results of streamflow change attribution are quite consistent. In fact, except for the Budyko-based methods denoted by Fu, Wang, and Yang, all other methods have shown a predominant role of human activity in the decline of streamflow in those catchments studied. More specifically, most methods revealed that the contribution of human activity to streamflow decline ranges from 71 to 79% with the corresponding contribution of climate ranging from 29 to 21%. By applying ER measures, there has been a greening trend across LP, especially in the last 20 years, and the slope gradient has been altered due to the water-conservation constructions [60,61]. The tradeoffs between conservation objectives and the negative impacts of ER measures on the streamflow remain to be addressed in the future works. Even though some methods (e.g., Fu, Wang, and Yang) are showing higher climate contributions compared to other methods (~50%), human activities still outweigh climate change in most catchment in terms of their impacts on streamflow, the results of which are consistent with previous contribution analyses in the same region [13,23]. The three methods (i.e., FU, Wang, and Yang), with which higher climate contributions would be derived, can be termed as parametric methods ( Table 2). They will typically put extra weights on factors such as topography, soil and vegetation conditions, and hydrological features at the catchment scale. However, considering more parameters may sometimes add more uncertainties, let alone geographical factors (e.g., location) and large-scale events, which can also affect the streamflow change. To draw conclusions on the sensitivity of various methods to parameters, elaborate research on a much larger sample size of catchments and corresponding hydrometeorological data is needed, but that is out of the research scope of this study. When comparing the performance of different Budyko-based methods at the catchment scale, the eight methods agreed to the largest degree in Jialu catchment with a standard deviation of 0.81%. On the other hand, the standard deviation of results derived by different methods in the Beiluo catchment reached 11.24%, which also indicates that further investigation is needed to determine the most suitable method to use when assessing human or climate contribution to streamflow change. The bar plots of standard deviation in Figure 3 reveal that the catchments can be divided into two categories based on the magnitude of the standard deviation. The first six catchments (i.e., Xinshui, Beiluo, Yanhe, Qingjian, Jinghe, and Dali) in Figure 3 have higher standard deviations, on average, compared to the rest of the catchments. By projecting the locations of catchments with the corresponding standard deviations of human/climate contribution (Figure 4), we can tell that those catchments are typically located in southern parts of LP and lower stream regions, compared to the rest of the catchments that are studied in this research. An exception is presented by methods using Hargreaves (HARG) E 0 estimation ( Figure A4), which corresponds to consistent low standard deviations across the catchments. Moreover, the results calculated by methods using the index-adjusted version of Hargreaves (HARGi) estimation show the highest standard deviations in the six austral catchments. More detailed discussions about the variability due to different E 0 estimations will be given in later sections.

Discrepancies Raised by Various E 0 Estimations
In Figure 5, for each of the 13 studied E 0 estimation methods, results were averaged over calculations using different Budyko-based methods. The same scenario of using color ramp to illustrate human and climate contributions to streamflow changes was applied as Figure 3 (i.e., the cells marked by dark red indicates that human factors predominantly contribute to streamflow change and the human contribution decreases as the cells become bluish). Figure 5 has confirmed the findings in Figure 3, that human activities are the primary driver for the streamflow decline in all the catchments studied in this research, although the magnitude of human contribution might vary due to the method applied and catchment characteristics. For the quantitative attributions using E 0 estimation methods, such as MAK, PT, and their coefficient-adjusted versions, the difference identified by the catchment location showed that boreal catchments (e.g., Kuye and Gushan) have higher human contributions to streamflow decline (~81%) than do austral catchments (e.g., Dali) (~74%). By calculating the standard deviation of contribution results at the catchment scale, different E 0 estimations tend to cause the highest variability in Xinshui, while it is the catchment Beiluo that has the highest variability of quantitative contribution in the case of comparing various Budyko-based methods (Section 3.2.1). Again, Jialu has witnessed the lowest variability in contribution results; the catchment Jialu may not be sensitive to either different Budyko-based methods or different E 0 estimations involved in this calculation. The quantitative attributions using E 0 , estimated by methods including FAO, JENSEN, and JENSENi, typically show reduced human contributions to the streamflow decline, while the contributions are still maintained at a level of 69%. The impacts of adjusting coefficients vary with different E 0 estimation methods; they either increase or decrease the final contributions calculated. However, for any particular coefficient adjustments, the sign of change (i.e., positive or negative) in the magnitude of contribution is always consistent in all catchments. For example, after adjusting the coefficient for the E 0 estimation method, ABTEW, the adjusted version (i.e., ABTEWi) leads to a decrease in human contribution to streamflow decline for all catchments. In other words, coefficient adjustments are independent from catchment characteristics such as geographical location. The bar plots of standard deviation in Figure 4 confirm the relatively lower sensitivity of the catchment Jialu to the method and E 0 estimation used in decomposing human and climate contributions. At least for austral catchments (e.g., Xinshui and Beiluo), the variabilities of contributions calculated by methods such as FU, WANG, and YANG are remarkably higher than those calculated by the other methods ( Figure A5). For catchments such as Xinshui, the standard deviation of quantitative contributions derived with different E 0 can be over 14% when using FU and WANG methods for calculation. At the catchment scale, the level of sensitivity of Budyko-based method to E 0 estimation also varies. E 0 has played a key role in calculating the sensitivity and contributions of climate change (which can also be partitioned into the contributions of P and E 0 separately) and human activities on streamflow variability with Budyko-based methods [1,13,21,23]. In the real world, the robustness of E 0 estimations by a particular method typically depends on the research domain (wet or dry climate) over which the method is applied. There is no doubt that a precise E 0 estimation would help to increase the accuracy of contribution analysis. Although identifying the optimal method to calculate E 0 for catchments on LP is not the focus of this study, the finding that, when using some estimation methods, the standard deviation is relatively lower for the contributions to the streamflow declines could provide implications for making further hypotheses regarding which method works the best on the selected catchments.

Uncertainties in the Quantitative Attribution
Typical sources of the uncertainty of quantitative contribution analysis include assumptions violating water balance assumption, the non-linearity of the interactions between climate change and human activities, as well as various ways to determine the value of E 0 . In terms of E 0 estimation, FAO has widely been used, while other methods could also be helpful when climate data does not feed the FAO's appetite. When examining the discrepancy shown by Figure 5, the HARG method formed the celling of human contribution of seven catchments, and after adjusting the coefficient, HARGi formed the floor of human contribution in the same seven catchments. Moreover, in the catchments Tuwei and Gushan, the situation reversed (i.e., HARG estimation led to the lowest human contribution, and HARGi estimation led to the highest estimation). Those results indicate that the different E 0 estimations that were involved in this research extended the uncertainty interval of quantitative contribution, which could be attributed to the deviation between the actual value and the estimation of E 0 . This study indicates that the uncertainty induced by E 0 estimation could weight a great share in total uncertainties. It is worth launching further sensitivity analyses upon the uncertainties in subprocesses of E 0 estimation, such as pre-process raw data (in situ ore remotely sensed data), equations for calculating E 0 , and interpolating point/station data, if necessary.
Surface runoff is an important part of the hydrologic cycle. This research attributed streamflow changes to the two general terms, climate change and human activities. However, it is the particular processes, such as changes in CO 2 concentration and land use and land cover (LULC), that directly or indirectly influence streamflow via biogeophysical and biogeochemical mechanisms. As discussed in Section 3.2.1, different Budyko-based methods greatly vary over the six catchments located at the austral part of LP, which may indicate that catchments with wet and warm climates are more easily subject to bias in contribution analysis of streamflow changes. However, a similar spatial pattern of contribution variability is not obvious for different E 0 estimations. This can also be demonstrated by mapping the standard deviation of quantitative contributions derived by a selected method (Budyko) and various E 0 estimation (Figure 6a) and the standard deviation of quantitative contributions derived by a selected E 0 estimation (FAO) with various calculation methods (Figure 6b). To address those uncertainties, LULC variables, such as normalized difference vegetation index (NDVI) and leaf area index (LAI), as well as meteorological variables, including temperature and precipitation anomalies, need to be investigated carefully.

Conclusions
In this research, changes in streamflow, evapotranspiration (E 0 ), and precipitation were first investigated for 12 catchments within Loess Plateau (LP) during the period of 1961-2018. With analysis methods such as the Mann-Kendall (M-K) test and the Range of Variability Approach (RVA), significant declining trends were found for streamflow in most catchments: low-flow seasons for the catchments studied were characterized by a decrease in E 0 ; and no significant trend was found for the change of precipitation in the study area on LP. Then, eight Budyko-based climate elasticity methods with 13 E 0 estimations were used to calculate the quantitative contributions of climate change and human activity to the decline of streamflow happening in those 12 catchments. With the whole period divided into baseline (the first 10 years) and evaluation (the remaining 48 years) periods, the results showed that human activity was the major contributor to the streamflow decline from baseline period to evaluation period for all 12 catchments, while the ratio of climate contribution to human contribution varied among catchments.
The discrepancies in results due to different methods and E 0 estimations applied were also examined. In conclusion, the sensitivity of a specific catchment to attribution methods or the key parameters in applied calculations vary and may be affected by geographical locations, landscape characteristics, and local climates. We argue that it would be difficult to find a combination of climate elasticity and E 0 estimation that is suitable for every case of contribution analysis, although it might be possible for catchments with similar climate types and catchment characteristics by more detailed sensitivity assessments. However, comparing common features and differences between methods as well as parameter values involved helps to capture the uncertainty in contribution analysis and could have implications for improving methods.

Acknowledgments:
The authors wish to thank the anonymous referees for their comments and suggestions. All individuals included in this section have consented to the acknowledgement.

Conflicts of Interest:
The authors declare no conflict of interest. Appendix A Figure A1. Time series of streamflow change for 12 catchments at intra-and inter-annual scales.  For Figures A1-A3, the time series for each time period is marked by different colors (blue lines for baseline period and yellow lines for evaluation period) and is fitted by linear regression model separately (subscript b refers to baseline period and subscript e refers to evaluation period).