Using Budyko-Type Equations for Separating the Impacts of Climate and Vegetation Change on Runoff in the Source Area of the Yellow River

Assessing the contribution rates of climate change and human activities to the runoff change in the source area of the Yellow River can provide support for water management in the Yellow River Basin. This paper firstly uses a multiple linear regression method to evaluate the contribution rates of climate change and human activities to the vegetation change in the source area of the Yellow River. Next, the paper uses the Budyko hypothesis method to calculate the contribution rates of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change of the Tangnaihai Hydrometric Station. The results showed that: (1) the annual runoff and precipitation in the source area of the Yellow River have a downward trend, while the annual potential evaporation and NDVI (Normalized Difference Vegetation Index) show an increasing trend; (2) The contribution rates of climate change and human activities to the vegetation change in the source area of the Yellow River is 62.79% and 37.21%, respectively; (3) The runoff change became more and more sensitive to changes in climate and underlying surface characteristic parameters; (4) The contribution rates of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change at Tangnaihai Hydrological Station are 75.33% and 24.67%, respectively; (5) The impact of precipitation on runoff reduction is more substantial than that of potential evaporation.


Introduction
In recent years, research on the mechanism of the hydrological cycle and water resources evolution under the changing environment has gradually become a hotspot [1]. Climate change and human activities are two important factors that lead to environmental change. Climate change directly affects the hydrological and thermal conditions of rivers and has a direct impact on water resources [2,3] and human activities have an indirect impact on the hydrological process of major basins, resulting in a series of water resources problems [4]. In recent decades, almost all river basins in the world have been affected by human activities and climate change at the same time. Vórósmarty [5] has shown that water resources have declined globally due to climate change and population growth and that this decline will continue.
The source area of the Yellow River provides more than 42% of the runoff of the entire basin [6] and is its main water production area and water conservation area. In recent decades, the runoff in the source area of the Yellow River has shown a gradual decrease, which has caused serious harm to the ecological security and economic development of the Yellow River basin. Therefore, discussing the characteristics of runoff change in the source area of the Yellow River and quantitatively assessing the impact mechanisms of climate change and human activities on the runoff change in the source area of the Yellow River can provide scientific insights and supports for the utilization of water resources in the Yellow River basin.
The Budyko hypothesis method is an effective method for statistically decomposing the changes in river runoff caused by climate change and human activities and has been widely used [7][8][9][10][11][12][13][14][15][16]. In recent years, some researchers have begun to use the Budyko hypothesis method to estimate the impact of climate change and human activities on the changes in runoff in the source area of the Yellow River. Zhou et al. [17] used the Budyko formula to quantify the impact of climate change and human activities on runoff change in the source area of the Yellow River during 1960-2010 and found that 1989 was the abrupt change year of runoff in the source area of the Yellow River basin and human activities are the main factor leading to runoff change, with a contribution rate of 67.62%. Zheng et al. [18] used the Budyko formula to determine the contributions of climate and catchment changes to runoff variation for every five years during 1961-2010 and found that climate change is the main factor leading to the runoff variation in the source area of the Yellow River, accounting for 60-70%. Wang et al. [19] used the Budyko formula to quantify the impact of climate change and human activities on runoff change in the source area of the Yellow River during 1960-2015 and found that the contribution of runoff reduction was 24.0% for precipitation change and 18.0% for potential evaporation. Further, for the Yellow River source region, the contribution of climate change scenarios and human activity factors to runoff change is close to 4:6. Shang et al. [20] used the Budyko hypothesis to quantify the impact of climate change and underlying surface change on runoff for every five years during 1956-2012 and found that the contribution rate of climate change to the runoff change of the Yellow River source was 52.47-68.50%.
In the Budyko hypothesis, runoff is mainly affected by precipitation (P), evaporation (ET), and underlying surface parameters (ω). The value of the underlying surface characteristic parameters (ω) depends on the soil type, topographical factors, and land types. Soil types and topographic factors generally do not change, so the values of underlying surface parameters are mainly affected by land type changes. However, previous studies often attributed land type changes to human activities. As a matter of fact, land type changes are affected by both climatic change and human activities. Therefore, existing researches often underestimated the impact of climate change on runoff change.
The grassland and forestland area in the source area of the Yellow River accounts for more than 80%, with building land and farmland area accounting for less than 1% [21]. Assuming that the soil type and topographical factors do not change in the source area of the Yellow River, we can use vegetation changes instead of land cover changes. NDVI (Normalized Difference Vegetation Index) data can effectively represent vegetation changes, so we used the NDVI data to characterize land type changes. To accurately assess the contribution rates of climate change and human activities to runoff change in the source area of the Yellow River basin, it is necessary to distinguish the impact of natural and human factors on land type changes, which is characterized by NDVI data.
Based on the above analysis, this paper first uses both the Mann-Kendall method and the cumulative anomaly method to determine the abrupt change year of river runoff over the source area of the Yellow River from 1961 to 2015 and next, uses a multiple linear regression method to calculate the contribution rates of climate change and human activities to the NDVI in the source area of the Yellow River. Finally, this paper uses the Budyko hypothesis method to calculate the contribution rates of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change in the source area of the Yellow River.

Study Area
The study area is the source area of the Yellow River, which is located in the upstream of the Tangnaihai hydrological station, between 95 • 50 -103 • 30 E longitude and 32 • 10 -36 • 05 N latitude. The elevation of the study area is higher in the west and lower in the east. The area with an elevation above 3000 m is 12.2 × 10 4 km 2 , accounting for about 17% of the total area of the whole Yellow River Basin. The reach is the main water production area and water conservation area of the whole basin. The annual average runoff of the study area is 232 × 10 8 m 3 and provides more than 42% of the total water resources in the Yellow River basin [6]. The study area is in the central Qinghai-Tibet Plateau. The daily temperature range is large, while the annual range is narrow. The annual average temperature ranges from −3.3 to 2.5 • C and the annual precipitation is between 332.4 and 743.2 mm. The geographic location of the source area of the Yellow River and the spatial distribution of meteorological stations are shown in Figure 1.

Study Area
The study area is the source area of the Yellow River, which is located in the upstream of the Tangnaihai hydrological station, between 95°50′-103°30′ E longitude and 32°10′-36°05′ N latitude. The elevation of the study area is higher in the west and lower in the east. The area with an elevation above 3000 m is 12.2 × 10 4 km 2 , accounting for about 17% of the total area of the whole Yellow River Basin. The reach is the main water production area and water conservation area of the whole basin. The annual average runoff of the study area is 232 × 10 8 m 3 and provides more than 42% of the total water resources in the Yellow River basin [6]. The study area is in the central Qinghai-Tibet Plateau. The daily temperature range is large, while the annual range is narrow. The annual average temperature ranges from −3.3 to 2.5 °C and the annual precipitation is between 332.4 and 743.2 mm. The geographic location of the source area of the Yellow River and the spatial distribution of meteorological stations are shown in Figure 1.

Data Sources
The data used in this study consist of three parts: (1) the observed annual runoff data of the Tangnaihai station from 1961 to 2015 are from the Yellow River Water Conservancy Commission (http://www.yrcc.gov.cn/); (2) The 20 meteorological station data in the source area of the Yellow River are obtained from China Meteorological Administration (http://www.cma.gov.cn). The dataset contains daily meteorological data of 20 stations in and around the study area from 1961 to 2015; (3) The NDVI data used in this study were obtained from the latest GIMMS NDVI3g version1 dataset released by the National Oceanic and Atmospheric Administration in 2016 (https://ecocast.arc.nasa.gov/data/pub/gimms/3g.v1/). Its time range is from July 1981 to December 2015 and its spatial resolution is 0.0833° × 0.0833°.

Data Sources
The data used in this study consist of three parts: (1) the observed annual runoff data of the Tangnaihai station from 1961 to 2015 are from the Yellow River Water Conservancy Commission (http://www.yrcc.gov.cn/); (2) The 20 meteorological station data in the source area of the Yellow River are obtained from China Meteorological Administration (http://www.cma.gov.cn). The dataset contains daily meteorological data of 20 stations in and around the study area from 1961 to 2015; (3) The NDVI data used in this study were obtained from the latest GIMMS NDVI3g version1 dataset released by the National Oceanic and Atmospheric Administration in 2016 (https://ecocast.arc.nasa.gov/data/pub/ gimms/3g.v1/). Its time range is from July 1981 to December 2015 and its spatial resolution is 0.0833 • × 0.0833 • .

Trend Analysis (Slope)
A simple linear regression model is the most commonly used method to analyze the trend change of observation data and has been widely used [22][23][24][25].
where n represents the total number of observation data; X i represents the value of the observation data at i year.

Mutation Analysis
A Mann-Kendall mutation analysis method and cumulative anomaly mutation analysis method are two widely used mutation test methods [26][27][28][29][30][31][32][33]. In order to identify the abrupt change year of runoff in the source area of the Yellow River as accurately and comprehensively as possible, these two methods were used to analyze the mutation year of runoff in this paper.
Mann-Kendall mutation test: A rank sequence is constructed for the time series matrix X containing n samples and the formulas are as follows: where UF 1 = 0, E(s k ) and Var(S k ) represent the mean and variance of S k , respectively. Their calculation formulas are as follows: Then, use the reverse order of X to perform the operation again according to the above formula (2)-(6). At the same time, make UB k = UF k (k = n, n − 1, . . . , 1), UB 1 = 0 to obtain the UB statistical curve. If the UF curve intersects with the UB curve and the intersection point is in the range of the 0.05 significance level, it can be considered that the year represented by the intersection point is the abrupt change year.
Cumulative anomaly method: for the time series matrix X containing n samples, the calculation formulas of cumulative anomaly values at time t are as follows: where X is the average value.

of 15
Then, the extreme point on the curve between the value of cumulative anomaly value S t and the time t is the mutation year.

Budyko Hypothesis
The water balance equation is expressed as [34]: where R represents runoff, P is precipitation, ET is evaporation, and ∆S represents the change of catchment storage water, which is considered 0 for a long period of time [35][36][37][38]. ET can be calculated based on the Budyko hypothesis. The formula is as follows [39,40]: where ET 0 represents the annual average potential evaporation (mm) and ω is the underlying surface characteristic parameters. ET 0 can be calculated using the Penman-Monteith equation.
where ∆ is the slope of the saturation vapor pressure versus temperature curve (kPa· • C −1 ), R n is the net radiation flux density at the surface (MJ/m 2 ·d), G is the sensible heat flux density from the surface to the soil (MJ/m 2 ·d), γ is the psychrometric constant (kPa· • C −1 ), U 2 is the wind speed at 2 m (m/s), T is the average temperature ( • C), e a is the saturation vapor pressure at air temperature (kPa), and e d is the actual vapor pressure of the air (kPa). Combining formulas (9) and (10), the water balance equation can be transformed into formula (12): The elastic coefficients of precipitation, potential evaporation, and underlying surface characteristic parameter to runoff can be calculated by formula (13)- (15) [15,34]: where ε P is the elastic coefficient of precipitation, ε ET 0 is the elastic coefficient of potential evaporation, and ε ω is the elastic coefficient of the underlying surface characteristic parameters.
The runoff data are divided into two periods: the base period (T 1 ) and the change period (T 2 ) based on the abrupt year of runoff at Tangnaihai hydrological Station in the Yellow River Basin. The annual average precipitation in the base period (T 1 ) is P 1 . The average annual precipitation in the changing period (T 2 ) is P 2 . The change in annual precipitation can be expressed as: The changes of potential evaporation (∆ET 0 ) and underlying surface characteristic parameter (∆ω) can be calculated from the T 1 period to the T 2 period using formulas (17) and (18): where ∆R P is the runoff variation value caused by annual average precipitation change from the T 1 period to the T 2 period, ∆R ET 0 is the runoff variation value caused by annual average potential evaporation change from the T 1 period to the T 2 period, and ∆R ω is the runoff variation value caused by underlying surface characteristic parameter change from the T 1 period to the T 2 period. The value of the underlying surface characteristic parameters (ω) depends on the soil type, topographical factors, and vegetation coverage. Assuming that the soil type and topographical factors do not change in the source area of the Yellow River, we can use vegetation changes instead of land cover changes. This paper uses the NDVI data to characterize vegetation changes and uses a multiple linear regression method to estimate the contribution rates of climatic change and human activities to vegetation changes [41]. First, the NDVI time series data are also divided into two parts based on the abrupt year of runoff data at Tangnaihai hydrological station: the base period (T 1 ) and the change period (T 2 ). Next, the average NDVI change in the two periods can be calculated as: where ∆NDVI is the average NDVI change in the two periods and NDVI T 1 and NDVI T 2 are the average NDVI values of the two periods, respectively. The change of NDVI can be attributed to the influence of natural factors (climate change) and human activities: where ∆NDVI C is the change of NDVI caused by climate change and ∆NDVI H represents the change of NDVI caused by human activities, such as afforestation, returning farmland to forest, etc. Several previous studies found that NDVI is closely related to precipitation and potential evaporation [41][42][43]. The multiple linear regression relationship between NDVI, the precipitation, and potential evaporation in the base period (T 1 ) can be estimated as follow: where a and b are regression coefficients and c is the intercept. By extending the formula (24), the NDVI value under climate change in the T 2 period can be calculated as follows: Water 2020, 12, 3418 7 of 15 The contribution of human activities and climate change to NDVI changes can be calculated separately: where NDVI T 2 ,S is the average NDVI under climate change in the T 2 period. η NDVI H and η NDVI C are the contribution rates of human activities and climate change to vegetation changes, respectively. The relative contribution rates of climate change (η R C ) and human activities (η R H ) to runoff change can be estimated using the following formulas:  From Figure 3a, the annual average precipitation in the source area of the Yellow River basin from 1961 to 2015 shows a slightly decreasing trend. The average decrease per year is 0.17 mm. As Figure 3b shows, there is an insignificant increasing trend of the annual average potential evaporation for the period 1961-2015, with an average increase of 0.45 mm per year.

Abrupt Changes of Runoff
In order to accurately identify the years with an abrupt change of runoff in the source area of the Yellow River basin, this study used both the Mann-Kendall mutation test and cumulative anomaly method to do an abrupt change analysis (Figures 4 and 5).

Abrupt Changes of Runoff
In order to accurately identify the years with an abrupt change of runoff in the source area of the Yellow River basin, this study used both the Mann-Kendall mutation test and cumulative anomaly method to do an abrupt change analysis (Figures 4 and 5).      Figure 5 shows that the cumulative anomaly values of runoff in the source area of the Yellow River basin from 1961 to 2015 first increased and then decreased, of which the turning point is 1989. Based on the results of the two different methods, this study believes that 1989 was the abrupt change year of runoff in the source area of the Yellow River basin.

Assessment of the Contribution Rates of Climate Change and Human Activities to Runoff Change
This study used precipitation (P), potential evaporation (ET0), and the underlying surface characteristic parameters (ω) to estimate the elastic coefficient of climate and land cover to runoff change. The underlying surface characteristic parameter (ω) depends on soil type, topographic

Assessment of the Contribution Rates of Climate Change and Human Activities to Runoff Change
This study used precipitation (P), potential evaporation (ET 0 ), and the underlying surface characteristic parameters (ω) to estimate the elastic coefficient of climate and land cover to runoff change. The underlying surface characteristic parameter (ω) depends on soil type, topographic factors, and vegetation coverage, assuming that the soil types and topographical factors in the source area of the Yellow River basin have not changed since 1961, and used vegetation change instead of land cover change. Previous studies often attributed LUCC (land use/cover change) to human activities. However, LUCC, particularly vegetation change, are impacted by both climate change and human activities. Therefore, previous studies are likely to underestimate the effect of climate change on runoff [11,15]. To accurately assess the contribution rates of climate change and LUCC to runoff change in the source area of the Yellow River basin, it is necessary to distinguish the impact of natural and human factors on LUCC.
Based on the results of abrupt analysis in Section 4.2, we dividend the period of 1961-2015 into two sub periods: the T 1 period  and the T 2 period . Several previous studies found that NDVI is closely related to precipitation and potential evaporation and can be an indicator to assess the impact of natural and human factors on LUCC [41][42][43]. The specific steps are as follows: (1) There were NDVI remote sensing image data from June 1981 to December 1989 in T 1 period , so we established the multiple linear regression relationship between NDVI and climatic factors (precipitation and potential evaporation) from July 1981 to December 1989. (2) Due to a lack of NDVI values in the period from January 1961 to June 1981, we further assumed that the multiple linear regression relationship in the period from June 1981 to December 1989 was the same as the period from January 1961 to June 1981. NDVI values for the period from January 1961 to June 1981 were calculated based on precipitation and potential evaporation in the corresponding periods. (3) We calculated the mean value of measured NDVI (NDVI T 1 ) in T 1 period , the mean value of measured NDVI (NDVI T 2 ) in T 2 period (1990-2015), and the mean value of NDVI under climate change (NDVI T 2 ,S ) in the T 2 period (1990-2015) according to the multiple linear regression function. (4) We calculated the contribution rates of climate change and human activities to NDVI changes using formula (26)-(29) (Table 1). Results showed that climate change plays a major role in increasing vegetation coverage and its contribution rate can reach up to 62.79%.
To further evaluate the long-term impact of climate change and the underlying surface characteristic parameter on runoff in the source area of the Yellow River, we used formulas (13)- (15) to calculate the elastic coefficients of precipitation, potential evaporation, and the underlying surface characteristic parameter in different years (Figures 6-8). As shown in Figures 6-8, the absolute values of ε P , ε ET 0 , and ε ω increased during the period 1961-2015, which means the runoff change became more and more sensitive to changes in climate change and the underlying surface characteristic parameter.      Table 2 showed the eigenvalues of different climatic and hydrologic variables in different periods: T1 and T2. Next, we calculated the differences in annual average precipitation, potential evaporation, and the underlying surface characteristic parameters between the two different periods. Then, the runoff variation value of Tangnaihai station caused by precipitation, potential evaporation, and LUCC were calculated by formula (19)- (21). Runoff variation caused by LUCC can be separated into two parts: natural and human factors; their contribution ratios can be obtained from Table 1.  Table 2 showed the eigenvalues of different climatic and hydrologic variables in different periods: T 1 and T 2 . Next, we calculated the differences in annual average precipitation, potential evaporation, and the underlying surface characteristic parameters between the two different periods. Then, the runoff variation value of Tangnaihai station caused by precipitation, potential evaporation, and LUCC were calculated by formula (19)- (21). Runoff variation caused by LUCC can be separated into two parts: natural and human factors; their contribution ratios can be obtained from Table 1. Based on the above method, the relative contribution rates of climate change (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to runoff at Tangnaihai station were calculated (Table 3). The contribution rates of climate change (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change of Tangnaihai Hydrometric Station were 75.33% and 24.67%, respectively. In general, climate change has played a major role in reducing runoff and the impact of precipitation on runoff reduction is more significant than that of potential evaporation.

Discussions
In the Budyko hypothesis, runoff is mainly affected by precipitation (P), evaporation (ET), and underlying surface parameters (ω). The value of the underlying surface characteristic parameters (ω) depends on soil type, topographical factors, and land types. Soil types and topographic factors generally do not change, so the values of underlying surface parameters are mainly affected by land type changes. The grassland and forestland area in the source area of the Yellow River accounts for more than 70%, while building land and farmland area account for less than 1% [21], so we can use vegetation changes instead of land cover changes.
NDVI data can effectively represent vegetation changes, so we used the NDVI data to characterize land type changes. To accurately assess the contribution rates of climate change and human activities to runoff change in the source area of the Yellow River basin, it is necessary to distinguish the impact of natural and human factors on land type changes, which is characterized by NDVI data.
Compared with previous studies, we first used NDVI data to characterize the land type change based on the actual situation of the study area in this paper and then distinguished the impact of natural and human factors on land type changes. In order to isolate climate impact on vegetation growth from human activities, a multiple linear regression method was introduced. Finally, the contribution rates of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change at Tangnaihai Hydrological Station were calculated. Results showed that previous studies underestimated the contribution rate of climate change on runoff change.
Although the contribution rate of climate change and human activities to runoff change is calculated, there are still some uncertainties: (1) Uncertainties of meteorological station observation data The observation data of meteorological stations are missing in some observation times. In addition, the precipitation and potential evaporation data of the study area are interpolated, but the results obtained by the interpolation still have some deviations from the actual distribution of precipitation and potential evaporation.
(2) Uncertainties in the division of base period and change period In order to identify the abrupt change year of runoff in the source area of the Yellow River as accurately and comprehensively as possible, Mann-Kendall and cumulative anomaly mutation analysis methods were used to analyze the mutation year of runoff in this paper, but we cannot guarantee that all mutation years were detected.
(3) Uncertainties caused by using NDVI to characterize land type change Building land in the source area of the Yellow River, which is attributed to human activities, accounts for less than 1%, but there are some uncertainties in using NDVI to characterize land type change.

Conclusions
This study first used the Mann-Kendall mutation test and accumulated anomaly method to identify the abrupt change year of runoff in the source area of the Yellow River. Next, a multiple linear regression method was applied to evaluate the contribution rates of climate change and human activities to the vegetation change in the source area of the Yellow River. Finally, we used the Budyko hypothesis to calculate the impact of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and human activities (mainly refers to vegetation changes caused by human activities) to the runoff change in the source area of the Yellow River.
Results showed that: (1) The contribution rates of climate change and human activities to the vegetation change in the source area of the Yellow River is 62.79% and 37.21%, respectively; (2) The runoff change became more and more sensitive to changes in climate and underlying surface characteristic parameters; (3) The contribution rates of climatic factors (including precipitation, potential evaporation, and subsequent vegetation changes) and vegetation changes caused by human activities to the runoff change at Tangnaihai Hydrological Station are 75.33% and 24.67%, respectively; (4) The impact of precipitation on runoff reduction is more substantial than that of potential evaporation.