Assessing the Impacts of Climate Change and Land Use / Cover Change on Runo ﬀ Based on Improved Budyko Framework Models Considering Arbitrary Partition of the Impacts

: Various models based on Budyko framework, widely applied to quantify the impacts of climate change and land use / cover change (LUCC) on runo ﬀ , assumed a ﬁxed partition used to distinguish the impacts. Several articles have applied a weighting factor describing arbitrary partitions for developing a total di ﬀ erential Budyko (TDB) model and a complementary Budyko (CB) model. This study introduces the weighting factor into a decomposition Budyko (DB) model and applies these three models to analyze runo ﬀ variation due to the impacts in the upper-midstream Heihe River basin. The Pettitt test is ﬁrst applied to determine a change point of a time series expanded by the runo ﬀ coe ﬃ cient. The cause for the change point is analyzed. Transition matrix is adopted to investigate factors of LUCC. Results suggest the consistency of the CB, TDB, and present DB models in estimating runo ﬀ variation due to the impacts. The existing DB model excluding the weighting factor overestimates the impact of climate change on runo ﬀ and underestimates the LUCC impact as compared with the present DB model. With two extreme values of the weighting factor, runo ﬀ decrease induced by LUCC falls in the range of 65.20%–66.42% predicted by the CB model, 65.01%–66.57% by the TDB model, and 64.83%–66.85% by the present DB model. The transition matrixes indicate the major factors of LUCC are climate warming in the upstream of the study area and cropping in the midstream. Our work provides researchers with a better understanding of runo ﬀ variation due to climate change and LUCC.


Introduction
Spatiotemporal variation of runoff has been an important component in hydrological cycle [1]. Climate change and land use/cover change (LUCC) are two major impacts on runoff [2]. Variation in temperature due to climate change causes the redistribution of precipitation to evapotranspiration and runoff [3]. Extreme hydrological events such as drought and flood influence runoff and intensify global water cycle [4][5][6]. On the other hand, LUCC such as deforestation and cultivation affects the regional water cycle and runoff [7]. A large amount of groundwater pumping for irrigation reduces runoff and increases evapotranspiration [8]. These therefore lead us to question how to differentiate and quantify the individual impacts of climate change and LUCC on runoff on reginal or global scale [9]. Exploring runoff variation due to the impacts helps researchers understand complex hydrological processes [10,11]. Inspecting the trend of hydro-meteorological series is a critical task for exploring the relations A great deal of effort has been made on the developments of the TDB and CB models with the weighting factor [21,25,26]. What seems to be lacking, however, is to examine the consistency of Water 2020, 12, 1612 3 of 15 multiple Budyko framework models based on the weighting factor accounting for arbitrary partitions of both climate change and LUCC impacts. In addition, a widely applied decomposition Budyko (DB) model has not considered the weighting factor [3].This study develops an improved DB model with the weighting factor and investigates the consistency of the DB, CB, and TDB models in assessing the impacts on runoff in the upper-midstream Heihe River basin during the study period of 1961-2014. The Mann-Kendall (MK) test and Sen's slope are adopted for trend analysis. The Pettitt test is applied to determine a reasonable change point of runoff coefficient time series. The Mezentsev-Choudhury-Yang function [22,27,28] is used to define the Budyko curve. The effect of the weighting factor on model predictions is explored. In addition, the transition matrix is adopted to investigate factors of LUCC from 1980 to 2000. Our work provides implications for not only better understanding of runoff variation but also more reliability of Budyko framework models for attribution analysis.

Mann-Kendall Test
Inspecting the trend of hydro-meteorological series is a critical task for exploring the relations between runoff Q and climate factors including precipitation P, and potential evapotranspiration E 0 . Sen's slope [29] and the non-parametric Mann-Kendall (MK) test [30,31] are used to identify temporal trends of P, E 0 and Q. Assume that each dataset of Q, P, and E 0 , expressed as {x i , x i+1 , . . . , x n } with n being the total number, is independent and identically distributed. The test statistic S and the variance of S are defined as: Var(S) = n(n − 1)(2n + 5)/18 The standardized test statistic Z for the standard normal distribution is calculated by: The significance level is set to 10%. With the significance level, the null hypothesis of no trend was rejected if |Z| > 1.645.

Sen's Slope
The slope β for a time series {x i , x i+1 , . . . , x n } proposed by Sen [29] is expressed as: A positive value of β indicates an increase, while a negative value indicates a decrease.

Pettitt Test for Determination of Change Point
A change point, dividing the study period into baseline and variation periods, needs to be determined before assessing the impact of LUCC on runoff. The baseline period before the change Water 2020, 12, 1612 4 of 15 point assumes no LUCC impact. The non-parametric Pettitt test is used to detect a change point of a time series [32]. Qiu et al. [33] revealed the Pettitt test to a time series of P or Q for the Heihe river basin did not give a unique and significant change point. For determining a change point, this study therefore applies the Pettitt test to runoff coefficient defined as σ = Q/P accounting for the relation between P and Q.
Considering a time series {x i , x i+1 , . . . , x N } with N being the total number, the Pettitt test uses the Mann-Whitney statistic U t,N that verifies if two subseries {x 1 , . . . , x t } and {x t+1 , . . . , x N } result from the same population. The test statistic U t,N is defined as: A most significant change point will be determined when the value of |U t,N | is a maximum, i.e., K t,N = max( U t,N ). The significance level associated with K t,N is expressed as: When ρ < 5%, the null hypothesis is rejected.

Budyko Framework
The Budyko [34,35] framework is based on the physical principles proposed by Schreiber [36] and Ol'Dekop [37]. The annual average actual evapotranspiration is dominated by the balance between the water supply of the atmosphere and the demand of atmospheric evaporation. A number of articles have presented various functions based on the Budyko framework for describing the relation between the long-term average of E/P and the dryness index of E 0 /P. The commonly used Mezentsev-Choudhury-Yang function is expressed as: where ω is a parameter reflecting catchment characteristic. By applying Equation (9) and ignoring change in the long-term storage, the runoff Q can be expressed as: The ω can thus be obtained by Equation (10) with annual values of observed P, E 0 , and Q. The change in observed runoff ∆Q obs is the difference in the mean annual runoffs Q b for baseline period and Q v for variation period, written as: The change in estimated runoff ∆Q est contains two components associated with the impacts of climate change ∆Q clim and LUCC ∆Q lucc , expressed as: Water 2020, 12, 1612

of 15
The contribution rates of climate change η clim and LUCC η lucc can be respectively written as [38]: Three models of CB, TDB, and DB are discussed in the following subsections.

Complementary Budyko (CB) Model
The CB model starts from complementary relationship of the elasticity coefficients associated with precipitation and potential evapotranspiration, expressed as [24]: Equation (15) reduces to: The change of runoff Q in Equation (16) is written as [25]: The weighting factor α is used to represent different paths from Point A to Point D in Figure 1 [21]. For Path 1 with α = 0, Point A to B reflects a fixed ratio of E 0 /P and a varying ω, indicating E 0 and P are assumed constant in baseline period and the catchment characteristic parameter is the only variable accounting for the impact of LUCC on runoff from baseline period to variation period. The path from Point B to D indicates the only impact of climate change on runoff with a fixed catchment characteristic parameter of variation period (i.e., ω v ) with subscript v being variation period. For Path 2 with α = 1, the segment from Point A to C accounts for the only impact of climate change on runoff from baseline period to variation period with a fixed catchment characteristic parameter of baseline period (i.e., ω b ) with subscript b being baseline period. The vertical descent from Point C to D reflects the only impact of LUCC due to the change of the catchment characteristic parameter (from ω b to ω v ). For an arbitrary path between Paths 1 and 2, Equation (17) with the weighting factor α becomes: The first two terms on the right-hand side stand for the impact of climate change, and the rest terms for the impact of LUCC.

Total Differential Budyko (TDB) Model
The TDB model with a first-order approximation to Q = P − f (P, E 0 , ω)P can be written as [22]: Similar to the derivation of Equation (18), the change of runoff for arbitrary paths represented by weighting factor α can be expressed as [25]: The first two terms on the right-hand side of Equation (20) are identical to those of Equation (18). The TDB model gives error in prediction because of applying the first-order approximation. This will be discussed in Section 4.2.

The Present Decomposition Budyko (DB) Model
Existing DB model is a graphic model based on Budyko curve represented by Path 2 with α = 1 (i.e., A-C-D) in Figure 1 [3]. With Path 2, the impacts of climate change and LUCC on runoff can be written as: where Due to the uncertainty of path [25], Path 1 with α = 0 (i.e., A-B-D) is also a possible way to assess runoff variation due to the impacts, expressed as: where (24), the present DB model with the weighting factor α accounting for an arbitrary path between Paths 1 and 2 can be expressed as: Table 1 shows the expressions of ∆Q clim and ∆Q lucc for the CB, TDB, and present DB models. Table 1. The expressions of climate change impact ∆Q clim and LUCC impact ∆Q lucc for the three models.

Transition Matrix of Land Use
Transition matrix of land use accounts for change of land use during different periods. The types of land use are classified into: crop-lands, forest, grassland, water bodies, snow and ice, urban and built-up, and barren in this study. The transition matrix of land use is obtained by using the Spatial Analyst Tools of ArcGIS function. One can refer to the study by Liu et al. [39] for detailed description of transition matrix.

Study Area and Datasets
The Heihe River, the second largest inland river in the south Qilian Mountain of northwestern China, suffers from serious water scarcity. The river generally flows northwards towards Mongolia. The main river channel (98 • -101 • E, 38 • -42 • N) is 821 km long. The total catchment area is 14.31 × 10 5 km 2 . Two main hydrological stations, Yingluoxia and Zhengyixia, divide the Heihe River basin into three sub-basins (upper, middle, and downstream). The upper and middle Heihe River basins, which contain alpine ice-snow and permafrost, mountainous forest zones, and a plain oasis agriculture zone, are selected as the study area ( Figure 2). The water sources of the upstream are mainly precipitation and glacier melt water. The midstream basin is an important commodity grain-producing area in a representative piedmont valley plain oasis.
Daily meteorological data at five meteorological stations in the study area were downloaded from the China Meteorological Data Sharing Service System (http://data.cma.cn/). Measurement of monthly runoff at the Zhengyixia hydrological station was collected from the Gansu Provincial Hydrological Bureau. The hydro-meteorological datasets for the upper and middle reaches spanned the period of 1961-2014. The daily potential evapotranspiration was estimated using the Penman-Monteith equation as suggested by Food and Agriculture Organization of the United Nations (FAO) [40]. The data of five meteorological stations with a resolution of 1 km × 1 km were spatially averaged by the Inverse Distance Weighted (IDW) method. The daily precipitation and potential evapotranspiration and monthly runoff were aggregated to obtain the individual time series expanded by 54 annual totals. In addition to the meteorological data, remotely sensed land use maps of 1980 and 2000 with a resolution of 1 km × 1 km were provided by the Resource and Environment Data Cloud Platform of the Chinese Academy of Sciences (http://www.resdc.cn/) and adopted to analyze the change in land use and cover in the study area.
The Heihe River, the second largest inland river in the south Qilian Mountain of northwestern China, suffers from serious water scarcity. The river generally flows northwards towards Mongolia. The main river channel (98°-101° E, 38°-42° N) is 821 km long. The total catchment area is 14.31 × 10 5 km 2 . Two main hydrological stations, Yingluoxia and Zhengyixia, divide the Heihe River basin into three sub-basins (upper, middle, and downstream). The upper and middle Heihe River basins, which contain alpine ice-snow and permafrost, mountainous forest zones, and a plain oasis agriculture zone, are selected as the study area ( Figure 2). The water sources of the upstream are mainly precipitation and glacier melt water. The midstream basin is an important commodity grainproducing area in a representative piedmont valley plain oasis.
Daily meteorological data at five meteorological stations in the study area were downloaded from the China Meteorological Data Sharing Service System (http://data.cma.cn/). Measurement of monthly runoff at the Zhengyixia hydrological station was collected from the Gansu Provincial Hydrological Bureau. The hydro-meteorological datasets for the upper and middle reaches spanned the period of 1961-2014. The daily potential evapotranspiration was estimated using the Penman-Monteith equation as suggested by Food and Agriculture Organization of the United Nations (FAO) [40]. The data of five meteorological stations with a resolution of 1 km × 1 km were spatially averaged by the Inverse Distance Weighted (IDW) method. The daily precipitation and potential evapotranspiration and monthly runoff were aggregated to obtain the individual time series expanded by 54 annual totals. In addition to the meteorological data, remotely sensed land use maps of 1980 and 2000 with a resolution of 1 km × 1 km were provided by the Resource and Environment Data Cloud Platform of the Chinese Academy of Sciences (http://www.resdc.cn/) and adopted to analyze the change in land use and cover in the study area.

Trend Analysis of Hydro-Meteorological Series and Determination of Change Point
The trend analysis of hydro-meteorological series is conducted by applying the MK test and Sen's slope for better understanding runoff process during the study period of 1961-2014. The results of trend analysis are shown in Table 2. The mean annual runoff (i.e., 28.47 mm/yr) is much smaller than the mean annual precipitation (i.e., 251.09 mm/yr) or potential evapotranspiration (i.e., 946.49 mm/yr). The trends of all meteorological factors were increasing. In contrast, a slight downward trend is detected for Q. The results of β > 0 or Z > 0 for P and β < 0 or Z < 0 for Q indicate a negative correlation between P and Q, which contradicts the linear positive correlation in the natural state. This may be attributed to the fact that irrigation and farming increased evapotranspiration and reduced runoff in the midstream of the Heihe River basin, important agricultural area [3]. It can be reasoned that the runoff variation is therefore related to both climate change and LUCC. The change point divides the study period of 1961-2014 into a baseline period and a variation period according to the Pettitt test [32]. Figure 3 shows temporal distributions of U t,N from the Pettitt test used to annual time series of precipitation P, runoff Q and runoff coefficient σ = Q/P.

The Impacts of Climate Change and LUCC on Runoff Based on the Budyko Framework Models
For better understanding the behavior of runoff change, the variation period of 1985-2014 is divided into six sub-variation periods of five years. Figure 4 illustrates the straight paths from the baseline period of 1961-1984 to the six sub-variation periods based on the Budyko framework. Each path moves from the baseline period to the upper left, indicating a greater increase in P than that in E 0 ; therefore, the climate therein became moister. The upper component of each path represents runoff increase due to climate change, and the leftward component accounts for runoff decrease caused by LUCC.
Water 2020, xx, x FOR PEER REVIEW 10 of 16 runoff increase due to climate change, and the leftward component accounts for runoff decrease caused by LUCC.  Table 3 displays the impacts of climate change and LUCC on runoff (i.e., ∆ and ∆ ) as well as the contribution rates (i.e., and ) estimated by the CB, TDB, and the present DB models considering the causes of α = 0 for the lower boundary of Path 1, α = 1 for the upper boundary of Path 2, and α = 0.5 in between during the variation period of 1985-2014. A reasonable range between two estimates for the cases of α = 0 and α = 1 can be seen. The estimates for the case of α = 0. 5 [22] function with the period of 1961-2014. In addition, our models consider the weighting factor α, but their model does not. Define difference D in the value of either ∆ or ∆ between Paths 1 for α = 0 and Path 2 for α = 1. Although the three models give close ranges of ∆ for the variation period of 1985-2014 in Table 3, the present DB model predicts a greater D than the others for each sub-variation period in Figure 5a. This may come from the fact that the present DB model, a graphical model, relies on two variables of ratio / and using Equation (9) while the others depend on three variables of , , and using Equations (18) and (20). It is worth noting that the CB and TDB models give exact the same D in Figure 5a because of the identical expression in Table 1 for ∆ . The TDB model, however, gives inaccurate results of much greater D for ∆ than the CB model for the first four sub-variation periods in Figure 5b because the first-order approximation of the TDB model causes a residual runoff change ∆ , defined as an observed runoff change ∆ minus an estimated one ∆ . Figure 6 shows the relations of ∆ and ∆ in panel (a) as well as ∆ and ∆ in panel (b) for the six sub-variation periods estimated by the TDB model with α = 0 for the lower boundary, α = 1 for the upper boundary, and α = 0.5 in between. When α = 0, ∆ is of  Table 3 displays the impacts of climate change and LUCC on runoff (i.e., ∆Q clim and ∆Q lucc ) as well as the contribution rates (i.e., η clim and η lucc ) estimated by the CB, TDB, and the present DB models considering the causes of α = 0 for the lower boundary of Path 1, α = 1 for the upper boundary of Path 2, and α = 0.5 in between during the variation period of 1985-2014. A reasonable range between two estimates for the cases of α = 0 and α = 1 can be seen. The estimates for the case of α = 0.5 fall in their own ranges. Existing DB model, a special case of the present DB model with α = 1, overestimates η clim and underestimates η lucc as compared with those predicted by the present DB model with α = 0.5. Reasonable ranges of 33.15%-35.17% for η clim and 64.83%-66.85% for η lucc are obtained using the present DB model. It is worth noting that the three models give close values of η lucc for the case of α = 0.5, indicating LUCC is the main cause for runoff decrease in the ranges of 65.20%-66.42% predicted by the CB model, 65.01%-66.57% by the TDB model and 64.83%-66.85% by the present DB model. Qiu et al. [33] reported η lucc = 53% predicted by the sensitivity model based on the Budyko framework. Our three models predict greater values of η lucc than theirs. The difference may result from the following reasons. Their model applies the Zhang et al. [12] function with the period of 1964-2006 while our models are based on the Yang et al. [22] function with the period of 1961-2014. In addition, our models consider the weighting factor α, but their model does not.
Define difference D in the value of either ∆Q clim or ∆Q lucc between Paths 1 for α = 0 and Path 2 for α = 1. Although the three models give close ranges of ∆Q clim for the variation period of 1985-2014 in Table 3, the present DB model predicts a greater D than the others for each sub-variation period in Figure 5a. This may come from the fact that the present DB model, a graphical model, relies on two variables of ratio E 0 /P and ω using Equation (9) while the others depend on three variables of E 0 , P, and ω using Equations (18) and (20). It is worth noting that the CB and TDB models give exact the same D in Figure 5a because of the identical expression in Table 1 for ∆Q clim . The TDB model, however, gives inaccurate results of much greater D for ∆Q lucc than the CB model for the first four sub-variation periods in Figure 5b because the first-order approximation of the TDB model causes a residual runoff change ∆Q res , defined as an observed runoff change ∆Q obs minus an estimated one ∆Q est . Figure 6 shows the relations of ∆Q est and ∆Q obs in panel (a) as well as ∆Q res and ∆Q obs in panel (b) for the six sub-variation periods estimated by the TDB model with α = 0 for the lower boundary, α = 1 for the upper boundary, and α = 0.5 in between. When α = 0, ∆Q res is of significantly positive relation with ∆Q obs . The |∆Q est | is underestimated as compared with |∆Q obs |. When α = 1, there is a negative correlation between ∆Q res and ∆Q obs , leading to an overestimate of |∆Q est |. Our results for the cases of α = 0 and 1 accord with the findings of Zhou et al. [25] and Wang et al. [26] that |∆Q est | by the TDB model was significantly underestimated for α = 0 and overestimated for α = 1 compared with |∆Q obs |. When α = 0.5, ∆Q res approaches zero and can be ignored, indicating the TDB model gives close values of ∆Q clim , ∆Q lucc , η clim , and η lucc to the others in Table 3. Error in ∆Q est by the TDB model, in other words, becomes minor for α approaching 0.5. Table 3. The impacts of climate change and LUCC as well as contribution rates calculated by the complementary Budyko (CB), total differential Budyko (TDB) and present decomposition Budyko (DB) models for three cases of α = 0 for the lower boundary, α = 1 for the upper boundary, and α = 0.5 in between.  Water 2020, xx, x FOR PEER REVIEW 11 of 16 significantly positive relation with ∆ . The |∆ | is underestimated as compared with |∆ |.
When α = 1, there is a negative correlation between ∆ and ∆ , leading to an overestimate of |∆ |. Our results for the cases of α = 0 and 1 accord with the findings of Zhou et al. [25] and Wang et al. [26] that |∆ | by the TDB model was significantly underestimated for α = 0 and overestimated for α = 1 compared with |∆ |. When α = 0.5, ∆ approaches zero and can be ignored, indicating the TDB model gives close values of ∆ , ∆ , , and to the others in Table 3. Error in ∆ by the TDB model, in other words, becomes minor for α approaching 0.5. Table 3. The impacts of climate change and LUCC as well as contribution rates calculated by the complementary Budyko (CB), total differential Budyko (TDB) and present decomposition Budyko (DB) models for three cases of α = 0 for the lower boundary, α = 1 for the upper boundary, and α = 0.5 in between.      Figure 7 shows the contribution rates η clim and η lucc for six sub-variation periods predicted by the CB, TDB, and the present DB models with α = 0 for the lower boundary, α = 1 for the upper boundary, and α = 0.5 in between. The three models agree to η clim and η lucc with minor error for each sub-variation period. The magnitude of η lucc is the largest in the sub-variation period of 2000-2004 for each model and each α and then decreases after the that period. This may be because China has experienced a striking economic boom with reform and opening-up policy since 1980 [38]. Population growth and economic development increased water withdrawal and utilization, resulting in a decrease in runoff. Since then, η clim generally increased and exceeded 50%, indicating that climate change impact dominates runoff change in the last two sub-variation periods. In addition, Figure 7 displays dramatically temporal variability of η clim and η lucc during the entire period except the last sub-variation period, which is similar to the finding of Wang et al. [26].

Factors of LUCC
A variety of studies have considered the change in ω (i.e., ∆ω) related to human activities. Roderick and Farquhar [41], however, revealed that ∆ω is linked to rainfall intensity of climate factor. Jiang et al. [42] estimated ∆ω according to the factors of temperature, potential evapotranspiration and irrigated area and found errors in ∆ω caused by assuming ∆ω due to human activities. Since LUCC is the dominant cause of runoff change in the study area as concluded in Section 4.2, this section applies the transition matrix of land use types for analyzing factors of LUCC. Table 4 shows the transition matrixes for describing LUCC from 1980 to 2000 in the upstream and midstream of the Heihe River basin. Grassland and barren occupied most of the area. In the upstream, the type of snow and ice was mainly converted into grassland and barren. Hao and Zong [43] also revealed this snow-ice area shrunk by half and turned into grassland and barren because of climate warming since 1980s. The reduction of the snow-ice area due to rising temperature and the spatiotemporal change of precipitation significantly affected runoff. In the midstream, four types (i.e., forest, grassland, water bodies, and barren) are converted into cropland, because the midstream is an important grain-producing area. More than eighty percent of water withdrawal was used to irrigate farmland, which is regarded as the main human activities in the midstream [44]. The increase in cropland indicates cropping is the major human activity resulting in LUCC. This result is consistent with the studies [45,46]. The grassland area also increased from the areas of the other types as a result of the Grain for Green Program. As concluded, the transition matrixes suggest the main factors of LUCC are climate warming in upstream and cropping in midstream, reflecting the impacts of not only climate change but also human activities. This result accords with the finding of Luo et al. [47] that runoff variation in the upstream is caused by climate change and in the midstream by LUCC.

Crop-Lands Forest Grass-Land Water Bodies Snow and Ice Urban and Built-Up Barren
Upstream

Conclusions
To conclude, this study is preliminary research on the impacts of climate change and LUCC on runoff in the upper-midstream Heihe River basin, but its relevance to Budyko framework models can be seen. A major finding is the consistency of the CB, TDB, and present DB models in predicting runoff variation due to the impacts for the study area. The results indicate the weighting factor α is applicable in providing a flexible and reasonable partition of both impacts on runoff. It can be reasoned that the impacts simultaneously induce runoff variation in the study area. Existing DB model excluding α gives an overestimated climate change impact on runoff and underestimated LUCC impact in comparison with the present DB model. In addition, error due to first-order approximation used to develop the TDB model becomes minor when α = 0.5. Furthermore, the extreme values of α = 0 and 1 can be used to quantify uncertainty of the impacts in providing those ranges in Table 3.
Our findings are consistent with those of the articles discussed above. With respect to the weighting factor α, our findings confirm those of Zhou et al. [25] and Wang et al. [26] although there are some differences regarding other aspects of the studies. These results lend some credence to the hypothesis that estimated runoff change |∆Q est | is underestimated by the TDB model for α = 0 and overestimates it for α = 1 as compared with observed runoff change |∆Q obs |. The discrepancy between ∆Q est and ∆Q obs results from the first-order approximation but becomes minor for α = 0.5. In addition, our finding is similar to that of Wang et al. [26]. The result indicates dramatically temporal variability of the contribution rates of climate change and LUCC to runoff during 1985-2014 because China has encountered striking economic boom since 1980. Moreover, our study agrees with Luo et al. [47]. The transition matrixes reveal the major factors of LUCC are climate warming in the upstream of the Heihe River basin and cropping in midstream. This finding is indicative of the fact that the reasons of inducing runoff variation are due to both climate change and LUCC in the study area.
This study has demonstrated the weighting factor is needed and should be stressed in assessing the impacts of climate change and LUCC on runoff. It follows that the use of multiple Budyko framework models are also needed for achieving reliability of model predictions. Our work provides researchers with a better understanding of runoff variation due to the impacts for semi-arid regions. However, whether this will also apply to humid regions of the world cannot be determined on the basis of this study. Further research is therefore warranted in exploring runoff variation induced by the impacts in humid regions using Budyko framework models with the weighting factor. Funding: The work has been supported by the grants from the National Natural Science Foundation of China (51809080), the Fundamental Research Funds for the Central Universities (2018B00114), and the National Key Technologies Research and Development Program (2018YFC0407900). All of the data in the study is available upon request. The authors would like to thank the Associate Editor, Prof. Wei, and two anonymous reviewers for their constructive comments and suggestions that greatly improve the manuscript.