## Abstract

Quantifying the contributions of climatic variables to runoff variance is still a great challenge to water resource management. This study adopted an extended Budyko framework to investigate the effects of terrestrial water storage changes (Δ*S*) on runoff variance across the source region of the Yellow River, China, during the period of 2003–2014. A new decomposition framework based on the extended Budyko framework was proposed to effectively quantify the contributions of different climatic variables including precipitation, *PET* and Δ*S* to runoff variance. The results demonstrated that the extended Budyko framework showed a better performance in presenting the water and energy balance than the original Budyko framework, especially at fine time scales. Meanwhile, the variance in runoff estimated by the new decomposition framework was close to that of runoff observations, indicating that this framework can effectively capture the variation in runoff during 2003–2014. It was also found that precipitation was the most important factor that contributed to runoff changes, while *PET* made a slightly smaller contribution compared to precipitation. Notably, the results also emphasized the important effects of Δ*S* on runoff variance at fine time scales, which was useful to better understand the interactions between atmospheric and hydrological processes for regions.

## HIGHLIGHTS

A new decomposition framework was proposed to effectively quantify the contributions of different climatic variables to runoff variance.

Precipitation was the most important factor that contributed to runoff changes, while

*PET*made a slightly smaller contribution compared to precipitation.The important effects of Δ

*S*on runoff variance at fine time scales can not be neglected.

## INTRODUCTION

Runoff is an important source of natural freshwater resources, which has been regarded as the fundamental element in social and economic development (Tikhamarine *et al.* 2020). Understanding the main drivers of runoff changes and variability is essential to decision-making in regional freshwater resource planning and management (Berghuijs *et al.* 2017). During the past decades, impacts of climatic variability on runoff magnitude have attracted widespread attention from many hydrologists and water resource managers. It is reported that about 31% of 145 major rivers across the world have shown significant changes in the mean annual runoff during the past decades (Walling & Fang 2003; Zhai & Tao 2017). Therefore, a comprehensive and quantitative understanding of runoff variance and its main sources is totally necessary for the evaluation of hydrological models (Wang & Hejazi 2011), climatic variability impact assessment (Collins *et al.* 2013) and uncertainty quantification in runoff estimation (Bock *et al.* 2018).

Numerous attribution methodologies have been proposed to separate the individual impact of climatic variability on runoff including hydrological modeling (Darvini & Memmola 2020), hydrological sensitivity method (Zuo *et al.* 2014) and climate elasticity method (Dey & Mishra 2017). For example, Ma *et al.* (2010) conducted a quantitative assessment of the influence of climatic variability on runoff reduction in the Miyun Reservoir basin, in China, by jointly using a distributed hydrological model and a climate elasticity model. Li *et al.* (2020) partitioned the contributions of precipitation and glacier melt to the runoff in a headwater basin of the Tarim River using a glacio-hydrological model and found that both rainfall and glacier melt are the primary causes of increased runoff in the study area during the period of 1971–2010. Chiew *et al.* (2006) estimated the sensitivity of runoff to the mean annual precipitation across the world based on an elasticity indicator and concluded that annual changes in precipitation were amplified in runoff, that is, a 1% change in the mean annual precipitation resulted in a 1–3% change in the mean annual streamflow globally.

Aside from various attribution methodologies mentioned above, the water balance method in combination with the Budyko framework has been considered as another useful approach to assess the contributions of different climatic variables and basin characteristics to runoff changes. The Budyko framework has been widely used in previous studies because it can accurately reflect the balance between water availability and energy supply across the study regions (Liang *et al.* 2014; Chen *et al.* 2020; Yang *et al.* 2020). Zhou *et al.* (2016) proposed a new partition method to distinguish the climate and basin attributes on the mean annual runoff within the Budyko complementary relationship. Zhang *et al.* (2019) conducted a Budyko-based framework with the purpose of quantifying the impacts of aridity index and other factors on annual runoff in 355 catchments selected from the United States. Notwithstanding that many attempts have been made to reveal the mechanism about how climatic variables influence runoff, there still exist some limitations in these studies. For example, the effect of changes in terrestrial water storage, which is an important component of the hydrological cycles for regions (Xie *et al.* 2019b), on runoff variance, is unknown because it has long been neglected when applying the Budyko framework to assess the effects of climatic variability on runoff. In addition, a quantitative analysis of variations in runoff at fine time scales (i.e. season or month) is also rarely made.

The source region of the Yellow River (hereafter SRYR), supplying water for millions of people, is an important region in China. Identifying and analyzing the variations in runoff across the SRYR is useful not only for understanding the water and energy balance between atmospheric and land surface hydrologic processes but also for improving the water resource management across the entire Yellow River basin (Xu *et al.* 2018). In addition, the mechanism behind runoff changes in the SRYR has attracted intensive concern in the decade owing to its link to the ‘Three Rivers Source Region Reserve’ and ‘Grain for Green Project’. Several studies have shown that the SRYR has experienced a significant decrease in runoff in recent years (Wang *et al.* 2018; Chu *et al.* 2019), which will further lead to variations in the temporal and spatial distributions of the water resources in the whole Yellow River basin. However, it has long been challenging to accurately assess the contributions of various natural drivers to runoff variance for this region due to its harsh environment and complex topography. Therefore, it is necessary and meaningful to quantify the contributions of various climatic variables to runoff, which will not only provide valuable guidance for decision-makers in making appropriate policies for water resource assessment but can also help hydrologists better understand runoff response to climatic variability.

As mentioned above, enough knowledge of the temporal variability of runoff is imperative to a comprehensive understanding of hydrological processes under climatic variability. Although previous studies have consistently investigated the effects of climatic variability on runoff to some extent, they fail to accurately quantify the contributions of these factors to runoff variance (Zarghami *et al.* 2011; Li *et al.* 2013; Lei *et al.* 2014). The main objective of this study is, therefore, to investigate the applicability of the extended Budyko framework by incorporating the effects of terrestrial water storage changes (hereafter Δ*S*) on regional water and energy balance at different time scales. Meanwhile, a new variance decomposition framework is developed to quantify the main factors influencing the runoff for the SRYR at intra-annual and inter-annual time scales, respectively. To our knowledge, current studies about the impacts of Δ*S* and other factors on runoff variance at fine time scales are still rare. Therefore, this study will be beneficial for accurately quantifying the response of hydrological cycles to climatic variability, especially for regions with limited hydrological data.

The other parts of this paper are mainly structured as follows. In the ‘Study area’ section, the study area is briefly described. The ‘Data’ and ‘Methods’ sections present the data and methods used in this study, respectively. The ‘Results’ section contains all the results illustrating the water balances of basins based on the extended Budyko framework at different time scales. Meanwhile, the effects of different factors on variations in runoff are also included in this section. Finally, all the discussions and major conclusions are presented in the ‘Discussions’ and ‘Conclusions’ sections, respectively.

## STUDY AREA

The SRYR mainly refers to the region controlled by the Tangnaihai (hereafter TNH) hydrological station located in the mainstream of the Yellow River, which is the second-longest river in China with a length of over 5,400 km (shown in Figure 1). The SRYR has a drainage area of 1.22 × 10^{5} km^{2} and accounts for approximately 15.2% of the area of the entire Yellow River basin. This region is situated in the northeastern Tibetan Plateau, which has an average altitude of over 4,000 meters above sea level. Due to its extremely harsh environment and climatic conditions, the SRYR has long been a sparsely populated region, and therefore this region can be viewed as a relatively pristine area with few human activities (Zheng *et al.* 2009). The SRYR, sometimes also termed as ‘the cistern of the Yellow River’, occupies only approximately 15.2% of the total area of the Yellow River basin but yields 35% of the total annual discharge of the Yellow River (Zheng *et al.* 2007; Yuan *et al.* 2018; Si *et al.* 2019). Therefore, the SRYR plays a considerable role in the downstream water resource management and hydrological cycles through the entire Yellow River basin.

The climate of the SRYR belongs to the typical Qinghai–Tibet Plateau climate system, which is characterized by distinct wet and dry seasons (Xu *et al.* 2018). Mainly influenced by the southwest monsoon from the Bay of Bengal, more than 70% of the annual precipitation across the SRYR occurs from June to September (Liang *et al.* 2010). For the period of 2003–2014, the mean annual precipitation of the SRYR is 583 mm while the annual *PET* is 837 mm, with an aridity index (*PET*/*P*) of 1.44. In addition, the annual runoff depth of the SRYR at the TNH hydrological station ranges from 109 to 220 mm with a mean value of 156 mm. Runoff coefficients range from 0.21 to 0.33 with an average of 0.27 during the same period mentioned above. More information about the related hydrological station, meteorological stations, geographical locations and drainage distribution of the SRYR have been depicted in Figure 1.

## DATA

### Evapotranspiration

Evapotranspiration (*ET*) is one of the most important elements in hydrological cycles, which can reflect the true state of water balance and energy exchanges in hydrological systems at different scales from field to global. However, *ET* has been extremely hard to measure or quantify through direct observations, especially in some regions with limited observations (Xing *et al.* 2018; Xie *et al.* 2019a). To estimate the variations in *ET* during the study period, three different datasets are jointly selected to estimate regional evapotranspiration in this study, namely (1) *ET* from the Global Land Data Assimilation System with Noah Land Surface Version 2 (hereafter ET_GLDAS, https://earthdata.nasa.gov/); (2) *ET* from Moderate Resolution Imaging Spectroradiometer (hereafter ET_MODIS, http://www.ntsg.umt.edu/project/mod16); and (3) *ET* from the Global Land Evaporation Amsterdam Model (hereafter ET_GLEAM, https://www.gleam.eu). Although these three *ET* products were derived from different principles or methods, all of them have been widely validated and applied in different regions including the SRYR because of their high quality and reliability (Xue *et al.* 2013; Liu *et al.* 2016; Ahlstrom *et al.* 2017; Asoka *et al.* 2017; Guo *et al.* 2017; Chang *et al.* 2018; Khan *et al.* 2018; Wang *et al.* 2018; Yang *et al.* 2021). In the following sections, the ensemble mean of these three independent *ET* products is used to describe the variations in actual *ET* across the study region.

### Meteorological data and *in situ* runoff observations

Meteorological data, including daily observations of mean air temperature (*T*), sunshine duration (*s*), relative humidity (*RH*) and wind speed (*u*) for the period 2003–2014, are collected from China Meteorological Administration (CMA, http://data.cma.cn/). All the above data are jointly used to obtain daily time series of potential evapotranspiration (*PET*) across the study region. In addition, the Thiessen polygon method is adopted to estimate daily time series of precipitation (*P*) across the study region using observations from different meteorological stations (shown in Figure 1), which are further aggregated to monthly results.

From the Yellow River Conservancy Commission, daily time series of runoff (*R*) observations during 2003–2014 are obtained at the TNH gauging station, which has been selected as the outlet of the SRYR. Furthermore, daily values of runoff observations are aggregated to monthly results with the purpose of keeping consistent with other variables in the water balance equation.

### GRACE-derived data

Terrestrial water storage mainly refers to all the forms of water stored under the earth (i.e. soil moisture and groundwater) and above the earth (i.e. lakes, rivers, wetlands and reservoirs). Variations in terrestrial water storage can reflect the budget between water input and water output for regions. The successful launch of the GRACE satellite has provided a new insight into the variations in terrestrial water storage, especially in some regions with few *in situ* observations (Song *et al.* 2015; Zhang *et al.* 2020). Since the GRACE satellite can effectively detect the variations in terrestrial water storage, in this study, three different GRACE mascon products are also jointly used to estimate Δ*S* with the purpose of validating the accuracy and reliability of results that are estimated by the water balance method. Note that some monthly values of GRACE data are missing because of ‘battery management’ during the study period. Therefore, these missing values are linearly interpolated from the previous and following months of the corresponding missing data. More detailed descriptions of datasets that are used in this study can be found in Table 1.

Variables . | Spatial resolution . | Temporal resolution . | Temporal extent . | Data references . | |
---|---|---|---|---|---|

Evapotranspiration (ET) | GLDAS2.0 Noah | 0.25° | Monthly | 1948–2014 | Rodell et al. (2004a) |

GLEAM | 0.25° | Daily | 2003–2017 | Martens et al. (2017); Miralles et al. (2011) | |

MODIS16 | 0.5° | Monthly | 2000–2014 | Mu et al. (2007) | |

Terrestrial water storage changes (ΔS) | GRACE-CSR Mascon | 0.5° | Monthly | 2002–2017 | Save et al. (2016) |

GRACE-JPL Mascon | 0.5° | Monthly | 2002–2017 | Landerer & Swenson (2012); Swenson & Wahr (2006) | |

GRACE-GSFC Mascon | 0.5° | Monthly | 2002–2017 | Awange et al. (2011); Luthcke et al. (2013) | |

Runoff (R) | – | Daily | 2003–2014 | – | |

Precipitation (P) | – | Daily | 2003–2014 | ||

Air temperature (T) | – | Daily | 2003–2014 | ||

Wind speed (u) | – | Daily | 2003–2014 | ||

Sunshine duration (s) | – | Daily | 2003–2014 | ||

Relative humidity (RH) | – | Daily | 2003–2014 |

Variables . | Spatial resolution . | Temporal resolution . | Temporal extent . | Data references . | |
---|---|---|---|---|---|

Evapotranspiration (ET) | GLDAS2.0 Noah | 0.25° | Monthly | 1948–2014 | Rodell et al. (2004a) |

GLEAM | 0.25° | Daily | 2003–2017 | Martens et al. (2017); Miralles et al. (2011) | |

MODIS16 | 0.5° | Monthly | 2000–2014 | Mu et al. (2007) | |

Terrestrial water storage changes (ΔS) | GRACE-CSR Mascon | 0.5° | Monthly | 2002–2017 | Save et al. (2016) |

GRACE-JPL Mascon | 0.5° | Monthly | 2002–2017 | Landerer & Swenson (2012); Swenson & Wahr (2006) | |

GRACE-GSFC Mascon | 0.5° | Monthly | 2002–2017 | Awange et al. (2011); Luthcke et al. (2013) | |

Runoff (R) | – | Daily | 2003–2014 | – | |

Precipitation (P) | – | Daily | 2003–2014 | ||

Air temperature (T) | – | Daily | 2003–2014 | ||

Wind speed (u) | – | Daily | 2003–2014 | ||

Sunshine duration (s) | – | Daily | 2003–2014 | ||

Relative humidity (RH) | – | Daily | 2003–2014 |

## METHODS

### Penman–Monteith method

*et al.*2008; Irmak & Mutiibwa 2010; Mcjannet

*et al.*2013; Mallick

*et al.*2015). The Penman–Monteith method is shown as follows:where

*PET*is the potential evapotranspiration (mm day

^{−1}); is the slope of the saturated vapor pressure–temperature curve (kPa °C

^{−1});

*R*

_{n}is net radiation at the canopy surface (MJ m

^{−2}day

^{−1});

*G*is the soil heat flux (MJ m

^{−2}day

^{−1}); is a psychrometric constant (kPa °C

^{−1});

*T*is the mean air temperature at a height of 2 m (°C);

*u*is the wind speed at a height of 2 m (m s

^{−1});

*e*

_{s}is the saturated vapor pressure at a height of 2 m; and

*e*

_{a}is the actual vapor pressure at a height of 2 m. Furthermore, daily time series of

*PET*are aggregated to monthly

*PET*with the purpose of maintaining consistency with

*ET*and

*P*.

### Estimation and validation of terrestrial water storage change (Δ*S*)

*S*on runoff variance have not been comprehensively examined before, mainly due to the lack of accurate observations of terrestrial water storage. To detect the changes in terrestrial water storage over the SRYR, the basin-scale water balance method is applied to estimate the monthly Δ

*S*across the SRYR (Wan

*et al.*2015; Lv

*et al.*2017), which can be described as follows:where Δ

*S*is terrestrial water storage changes (mm);

*P*is precipitation (mm);

*R*is runoff (mm);

*ET*is evapotranspiration (mm), which can be estimated by the average of three different evapotranspiration products as mentioned in the ‘Evapotranspiration’ section.

*S*across the SRYR are obtained based on the water balance method (Equation (2)). To further validate the accuracy of these results derived from the water balance method, the changes in terrestrial water storage across the SRYR are also estimated based on GRACE data independently:where is GRACE-derived terrestrial water storage changes (mm);

*TWSA*(

*m*+ 1) is terrestrial water storage anomalies (mm) for month (

*m*+ 1), while

*TWSA*(

*m*− 1) is terrestrial water storage anomalies (mm) for month (

*m*− 1). Previous studies (Long

*et al.*2014; Xie

*et al.*2019a) have indicated that this method can not only be used effectively to estimate the regional Δ

*S*but also includes some light numerical smoothing.

### Original and extended Budyko frameworks

*ET*for regions can change under the joint effects of energy and water availabilities. Hence, the Budyko framework is reflecting water and energy balance for regions, which can be represented as Fu's equation (Fu 1981; Zhang

*et al.*2004):where is the long-term average aridity index for study regions;

*P*and

*PET*are precipitation (mm) and potential evapotranspiration (mm), respectively; is the parameter reflecting basin-specific characteristics, such as soil moisture, vegetation cover and climate seasonality (Yang

*et al.*2008; Zheng

*et al.*2018), which can be estimated by the least squares method.

*et al.*2013; Greve

*et al.*2020). In fact, recent studies have also made attempts to extend this equation with the purpose of validating the variability of annual or even monthly water balance in different regions (Yang

*et al.*2007; Zhang

*et al.*2008; Carmona

*et al.*2014). Therefore, Fu's equation at monthly or annual scales has been further modified as:where ∅︀ (=

*PET*/

*P*) and

*n*represent the aridity index and the basin characteristics parameter at different time scales such as annual or monthly.

*et al.*2009; Istanbulluoglu

*et al.*2012). Wang (2012) pointed out that the influences of Δ

*S*should not be neglected anymore when applying the Budyko framework at annual or monthly scales for regions. That is to say, Δ

*S*may play a more important role in determining hydrologic response at the intra-annual scales. To better describe the actual variability of water balance for regions under non-steady-state conditions, it is very necessary to take the influences of Δ

*S*into consideration, which can be obtained from Equation (2). Chen

*et al.*(2013) suggested replacing atmospheric water supply (

*P*) by the total available water (

*P′*) in Equation (5), and therefore the original Budyko framework can be extended into:where

*P′*(=

*P*− Δ

*S*) is effective precipitation (mm), which represents the total available water for regions and mainly depends on both atmospheric water supply and basin storage (Wu

*et al.*2017a,); ∅︀

*′*(=

*PET*/

*P′*) represents the aridity index considering Δ

*S*; and

*n*is the basin characteristics parameter, which can be obtained based on the least squares method.

### Runoff variance decomposition framework

*P*and

*PET*) and basin storage changes (i.e. Δ

*S*) on runoff, a new variance decomposition framework is proposed in this study. According to the water balance equation and assuming that the long-term storage change , the observed runoff (

*R*) deviation from its long-term mean () at a specific time interval can be expressed as follows:where , , and refer to the changes in runoff, precipitation, evapotranspiration and terrestrial water storage at a specific time interval, such as month, season or year.

*ET*(i.e. ) for regions mainly consist of different changes in

*P*,

*PET*and Δ

*S*, respectively, which can be described as:where , and are the same as that shown in Equation (7), represents changes in potential evapotranspiration, is the long-term average aridity index, and are the Budyko framework and its first-order derivative, respectively.

*R*can therefore be derived by taking the square of Equation (9), summing over

*N*samples and scaled by

*N*− 1:where and indicate the variance and covariance of different factors, respectively. Additionally,

*w*before the different variance (or covariance) terms represents the corresponding weighting factors, respectively, which can effectively quantify the contributions of different factors to

_{i}*R*variance and can be analytically estimated from the long-term average aridity index () in combination with the Budyko framework. Different weighting factors are presented as follows:where is the long-time average aridity index; and represent the original Budyko framework and its first-order derivative, respectively; and

*w*refers to different weighting factors. In these equations, a positive (or negative) weighting factor indicates that an increase (or decrease) in the corresponding term (such as , and ) will result in an increase (decrease) in runoff variance . For example, if

_{i}*w*is equal to 0.1, this means that a 10% increase of precipitation variance would simultaneously bring a 1% increase in runoff variance . According to Equations (10) and (11)–(16), runoff variance across the SRYR can be finally partitioned into variances from climatic variables and basin storage including

_{p}*P*,

*PET*and Δ

*S*.

## RESULTS

### Monthly time series of *P*, *PET*, *ET* and *R* during 2003–2014

Figure 2 presents the monthly time series of *P*, *PET*, *ET* and *R* during the study period across the selected SRYR. It has been challenging to acquire accurate *in situ* measurements of actual *ET* for large basins such as the SRYR with few meteorological and hydrological observations. As shown in Figure 2, the mean monthly *ET* across the SRYR ranges from 11.5 to 88.2 mm during 2003–2014. In addition, monthly *PET* shows distinct seasonal variations with a wide span from 23.5 to 123.8 mm. Mean monthly precipitation has a wide range from 0.3 to 188.5 mm. Generally, it shows an obvious seasonal cycle with a maximum value in summer (June to August) and a minimum value in winter (December to February), which is in line with the results from previous studies (Meng *et al.* 2014; Deng *et al.* 2020). Monthly time series of runoff is also acquired from the TNH station, which shows a reasonable correspondence with the precipitation and *ET* shown in Figure 2.

### Validation of changes in terrestrial water storage across the SRYR

To better understand the variations of water cycles and their response to climatic variability across the SRYR, the Δ*S* at different time scales (monthly, seasonal and annual, respectively) is needed. In this study, regional Δ*S* is estimated as the residual of water balance closure. As shown in Figure 3, monthly time series of Δ*S* based on the water balance method shows a seasonal variation with a wide range from −41.2 to 59.8 mm.

Since the aforementioned results are directly estimated by the water balance method, it is necessary to further validate monthly time series of Δ*S* against the GRACE-derived observations across the SRYR. The values of RMSE and correlation coefficient (*r*) are 14.05 mm and 0.79, respectively. Overall, there exists a statistically significant positive correlation at the 0.05 significance level (*p* < 0.05) between Δ*S* estimated by the water balance method and that derived from GRACE data across the study region, which reach to maximum and minimum values almost simultaneously during 2003–2014. All the results shown in Figure 3 indicate that monthly time series of Δ*S* based on the water balance method are reasonable and acceptable for the SRYR. Furthermore, estimated monthly Δ*S* can be added into annual and seasonal time series of Δ*S*, respectively.

### Effect of Δ*S* on water balance at different time scales

#### Monthly scales

The Budyko curves for the monthly datasets across the SRYR during 2003–2014 have been shown in Figure 4. *ET* and *PET* for regions are scaled by *P* and *P′*, respectively, with the goal of illustrating the importance and necessity of considering Δ*S* when applying the Budyko framework to assess runoff variance over a time interval. Due to the influence of Δ*S* on hydrological cycles, there exist significant differences between Figure 4(a) and 4(b). As shown in Figure 4(a), the evapotranspiration ratio (*ET*/*P*) increases linearly with the increasing aridity index (*PET*/*P*) when neglecting the effects of Δ*S* at the monthly scale. In addition, many data points of *PET*/*P* and *ET*/*P* obviously depart from the theoretical Budyko curve, and some scatter points even fall above the water limit line of Budyko space (i.e. *ET*/*P* >1). This is inconsistent with the previous results found by Chen *et al.* (2020), who stated that the obvious linear relationship between *ET*/*P* and *PET*/*P* can be captured at the monthly scale in some typical arid and semi-arid regions.

In contrast, the data points of *PET*/*P′* and *E*/*P′* considering Δ*S* generally follow the theoretical Budyko curve well when considering Δ*S* for the SRYR at the monthly scale as shown in Figure 4(b). Overall, the points scaled by *P′* (i.e. with considering Δ*S*, Figure 4(b)) are denser than those scaled by *P* (i.e. without considering Δ*S*, Figure 4(a)) at the monthly scale. Additionally, the data points considering Δ*S* (*R*^{2} = 0.49) show a better performance in the fitting of Budyko curves than those without considering Δ*S* (*R*^{2} = 0.11), indicating that Δ*S* plays an important role in the basin-scale water and energy balance. The above results further reflect that regional Δ*S* indeed makes up a large proportion of the partitioning of *P* into *R* (or *E*), and therefore it should be taken into consideration when analyzing the hydrological cycles at the monthly scale.

#### Seasonal scales

Investigating the water and energy budget over regions at seasonal scales can help us better understand the variability in runoff, which is important and meaningful for predicting some extreme events such as droughts or floods in advance under climatic variability, especially in some ungauged basins. Therefore, the Budyko framework is also applied in the SRYR using seasonal time series of data during the study period. Similar to the results shown in Figure 4, Figure 5 also presents the water balance in the Budyko curve for the SRYR. The results in Figure 5(a) demonstrate that the data points representing *PET*/*P* versus E/P distribute as linear curves rather than the Budyko curves as expected while assuming that Δ*S* was negligible at the seasonal scale. Some data points exceeded the ‘water limit’ boundary represented by *ET*/*P* = 1 in Figure 5(a), indicating that actual *ET* is more than the amount of *P* due to the neglect of Δ*S* in the water balance equation.

Figure 5(b) presents the ratio of seasonal actual evapotranspiration to available precipitation (*ET*/*P′*) as a function of corresponding dryness index (*PET*/*P′*) considering Δ*S* for the SRYR. As expected, the data points considering Δ*S* in the study period present a perfect Budyko relationship at the seasonal time scale as depicted in Figure 5(b). In other words, the extended Budyko framework considering Δ*S* (i.e. Equation (6)) can accurately capture the seasonal variability of the regional water and energy balance across the SRYR in this figure.

#### Annual scales

At the annual scale, both the scatter points shown in Figure 6(a) and 6(b) follow an approximate Budyko-like distribution as expected, with none of them falling in the region limited by water and energy boundary. Furthermore, it can be found that most of the data points are well located in the nearby region of the theoretical Budyko curve except one point shown in Figure 6(a), which may be caused by the neglect of Δ*S*.

Figure 6(b) shows the corresponding extended Budyko framework for the annual datasets across the study region during the study period. The results shown in Figure 6(b) show that all data points representing *PET*/*P′* versus *ET*/*P′* distribute in the regions of theoretical Budyko curve perfectly. Namely, the balance between water availability and energy supply at the annual scale for the SRYR can be explained well by the extended Budyko relationship. This finding is highly in line with Yang *et al.* (2007) who stated that this extended Fu's equation can be used for predicting the inter-annual variability of regional water balances.

In general, using the newly modified effective precipitation (*P′* = *P* – Δ*S*) as the proxy of total water availability for closed basins, the data points of *PET*/*P′* and *ET*/*P′* (*R*^{2} = 0.81) generally show a better performance in the fitting of the Budyko relationship than those without considering Δ*S* (*R*^{2} = 0.67). In particular, the original Budyko framework using the precipitation as the available water resource is reasonable and applicable for the long-term hydrological cycles of natural and closed catchments because Δ*S* is relatively small compared to variations in runoff and the other variables such as *P* or *ET* at the annual scale. When the time scale becomes finer, the contribution of Δ*S* to the water balance will become bigger, which indicates that the effects of Δ*S* are generally significant in the SRYR at the monthly and seasonal scale.

### Contributions of climate factors and basin storage to runoff variance based on the extended Budyko framework

It is widely known that runoff is closely related to regional water resource management and planning, which is of critical importance for sustainable social and economic development (Mccabe & Wolock 2011; Fowler *et al.* 2016). In fact, variation in runoff is usually viewed as one of the most important indicators that can reflect the true state of regional available water. Therefore, deep insights into runoff variance can help us better assess how regional water availability (in the form of runoff) has changed in the past decade. According to the newly proposed variance decomposition framework (i.e. Equation (10)), a quantitative assessment of the contributions of different factors in hydrological cycles (i.e. *P*, *PET* and Δ*S*) to the variance in runoff at different time scales can be implemented across the SRYR. By incorporating Δ*S* into the Budyko framework, the critical role of Δ*S* in water and energy cycle dynamics, especially in runoff variance, has been highlighted. Results of different weighting factors shown in the variance decomposition framework (i.e. Equation (10)) are shown in Figure 7, which represents the corresponding contributions to runoff variance when the same changes occurred in each term. It can be found that the value of *w _{P}* is identical to that of

*w*

_{ΔS}shown in Equation (10), which indicates that

*P*and Δ

*S*can make same contributions to runoff variance on the condition that same variations occur in

*P*and Δ

*S*, respectively. Additionally, the negative contribution from the covariance between

*P*and Δ

*S*(i.e.

*w*Δ

_{P,}*S*) to runoff variance is the most significant among all weighting factors with a value of −0.43. In comparison with the other weighting factors listed in Figure 7, the contribution of variance in

*PET*(i.e.

*w*) to runoff variance is relatively insignificant with a value of +0.03.

_{PET}The contributions of each variable to runoff variance at different time scales have been depicted in Figure 8. As shown in Figure 8(a), runoff variance at the monthly scale mainly comes from *P*, while the contribution from *PET* is relatively small. Notably, Δ*S* obviously plays an important role in runoff variance when compared to *PET*, although it has been neglected in previous studies (Jiang *et al.* 2015; Wu *et al.*2017b; Wang *et al.* 2019). Some other terms in Equation (10), such as and , generally make obviously negative contributions to runoff variance, which indicate that the source of runoff reduction mainly arises from the covariance between *P* and *PET* (or Δ*S*). In addition, variance in runoff estimated by Equation (10) has a good agreement with that derived from observations at the monthly scale with a small bias of 7%, which demonstrates that the variance decomposition framework is applicable and effective to assess the contributions of terms to runoff variance across the SRYR.

The contributions of different factors to variance in runoff at the seasonal scale show a similar pattern to that at the monthly scale after comparing Figure 8(a) and 8(b). As shown in Figure 8(b), *P* is still the most important source that results in the variations in runoff. At the same time, both Δ*S* and *PET* also have positive effects on runoff variance. In contrast, the negative contribution to the runoff variance results from the terms, including and , which is similar to that shown in Figure 8(a). Additionally, variance in runoff predicted by Equation (10) is 847 mm^{2}, which is rather close to the results derived from observations at the seasonal scale (i.e. 736 mm^{2}) with a relative bias of 15%.

The contributions of different climate factors and basin storage to runoff variance at the annual scale are also plotted in Figure 8(c). Undoubtedly, *P*, as a considerably remarkable impact factor, greatly affects runoff variance at the annual scale. Δ*S* is another critical contributing factor in enhancing the inter-annual variance in runoff. Notably, it can be observed from Figure 8(c) that the contribution of Δ*S* to runoff variance is more significant at the monthly scale than that at the annual scale because Δ*S* is more sensitive to the monthly climatic fluctuation than annual fluctuation. In addition, the terms of covariance between *P* and Δ*S* have a negative contribution to variance in runoff. The contributions of other terms to inter-annual runoff variance also have been shown in Figure 8(c), respectively.

According to the proposed variance decomposition framework, the contributions of different climate factors and basin storage to runoff variance for the SRYR have been fully quantified. Overall, *P* contributes more to runoff variance than the other factors such as *PET* and Δ*S* because it is the most important input for regional water balance. Especially, this study emphasizes the effect of Δ*S* on runoff variance. As mentioned in Figure 8, the effects of Δ*S* on runoff variance mainly come from three parts, that is, Δ*S* and the related covariance terms including and . It should be noted that Δ*S* can contribute more to the variance in runoff than *PET*. Although the effects of Δ*S* have long been neglected in previous studies (Jiang *et al.* 2015; Wang *et al.* 2019), the results from this study have demonstrated that Δ*S* may play a critical role in runoff variance. Obviously, the runoff variance simulated by Equation (10) would be possibly underestimated or overestimated at any time scales without considering the effect of Δ*S*.

## DISCUSSIONS

### Performance of Budyko framework at different time scales

To investigate the role of Δ*S* (including variations in surface water storage, soil moisture storage and groundwater storage) in the water balance of the SRYR, the variation in terrestrial water storage has been firstly estimated, which shows good agreement with that derived from GRACE data. Δ*S* and other factors including *P* and *PET* were jointly used to analyze the regional water and energy balance at different time scales based on the Budyko framework as mentioned in the ‘Results’ section.

Overall, the results demonstrate that it is reasonable to exclude the variations in terrestrial water storage from precipitation (i.e. *P′**=**P* − Δ*S*), which is more consistent with the water availability concept in boundary conditions of the original Budyko framework, especially at finer time scales such as monthly or seasonal (Chen *et al.* 2020). Meanwhile, this result repeatedly shows the importance of Δ*S* when analyzing the regional water balance, since the water available estimation by using *P* only is highly overestimated. Similar results also have been obtained in some other regions around the world (Zeng & Cai 2016).

### Main factors controlling runoff variance at different time scales

The variance decomposition framework of the runoff was applied at monthly, seasonal and annual scales independently. The results show that the runoff variance simulated by the decomposition framework is in line with that observed by the hydrological station in general. However, there still exists some discrepancy between the simulated runoff variance and that derived from observations, especially at the annual scale. This phenomenon can be explained by the following two factors. On one hand, precipitation is viewed as the sole water input for hydrological cycles in the study region when estimating Δ*S* based on the water balance method, which has proved to be reliable and robust across the SRYR. However, some other natural water resources, such as the melting of snow or glaciers, are not taken into consideration due to the lack of *in situ* observations. On the other hand, the term of aridity index (i.e. ) in Equations (11)–(16) is supposed to reflect the long-term mean climate condition for the SRYR while the annual sample data in this study are relatively limited, which is an important reason why the estimated runoff variance is considerably lower than the observations at the annual scale when applying the decomposition framework. The decomposition framework proposed in this study can be better validated on the condition that more observations about Δ*S* at the annual scale are available.

### Uncertainties and limitations

In this study, different terms in the water balance method are realistically derived from the meteorological stations. In fact, it remains a big challenge to obtain accurate and reliable estimations about precipitation or *PET* in high-altitude regions such as the SRYR. Furthermore, extremely limited meteorological stations also make it more difficult to capture the true state of water and energy balance cycles in this region. *ET*, as one of the most difficult variables to obtain or measure, at regional and basin scales, is difficult to be accurately estimated due to its link to Earth's water, energy and carbon cycles (Wang *et al.* 2018). Therefore, the average of three widely used *ET* products is adopted to describe the variations in *ET* with the goal of reducing error and uncertainty. However, the uncertainty of different terms mentioned above will further lead to an error in Δ*S* through the water balance method (i.e. Equation (2)) via the principle of uncertainty propagation. As documented in previous studies (Rodell *et al.* 2004b; Xie *et al.* 2019a), Δ*S* estimated by the water balance method is vulnerable to input uncertainty induced by different terms as mentioned above. Although the observations derived from the GRACE satellite can further help us validate the accuracy and reliability of Δ*S* estimated by the water balance method to some extent, it still inevitably leads to some errors or uncertainties in the estimation of Δ*S* for the SRYR due to its complicated geophysical conditions and extremely limited observations. In other words, the input uncertainty of Equation (10) may be an important reason why different data points (*ET*/*P* versus *PET*/*P* or *ET*/*P′* versus *PET*/*P′*) are not fully fitted with the Budyko curve.

The newly proposed variance decomposition framework is physically robust, while some high-order terms in this equation have been neglected in this study for simplicity. Although the errors induced by these high-order terms are relatively limited, it may still lead to some discrepancies between the estimated and theoretical results. Given the above reasons, runoff variance estimated by the proposed variance decomposition framework would be possibly overestimated or underestimated. Therefore, more efforts will be made in our next study to reduce these errors and obtain more reliable results.

## CONCLUSIONS

This study applied the original and extended Budyko framework to investigate the role of Δ*S* in water and energy balance across the SRYR at different time scales. A new variance decomposition framework was also proposed to partition the variance in runoff into different climate factors and basin storage. The major findings from this study are summarized as follows:

- (1)
Variations in terrestrial water storage have proved to play an important role in the hydrological cycles, and neglecting the effects of Δ

*S*would result in obvious errors when applying the original Budyko framework to analyze the water and energy balance for regions, especially at fine time scales such as month or season. In comparison with the original Budyko framework, the extended Budyko framework considering Δ*S*can better reflect the true process of exchange between water and energy across the study region. - (2)
The variance decomposition framework proposed in this study can be applied to effectively assess the contributions of different climatic factors and basin storage to the variance in runoff at different time scales. It was found that this method is more effective at monthly or seasonal scales than at the annual scale, which may result from the limitation of sample data available at annual scales.

- (3)
According to the proposed variance decomposition framework,

*P*is the main source of the variance in runoff for the SRYR and has significantly positive effects on variance in runoff. Δ*S*makes a positive contribution to variance in runoff, which is more significant at monthly and seasonal scales than the annual scale. In addition, the different covariance terms between*P*and*PET*(or Δ*S*) can have a significant influence on variance in runoff.

The conclusions drawn from this study may help us better understand the response of hydrologic cycles under climatic variability across the SRYR, which can provide valuable guidance for decision-makers in evaluating and predicting water resources. According to the extended Budyko framework, the variance decomposition framework proposed in this study can make a quantitative assessment of the influence of climatic variability on runoff changes, which can be applied in other study regions.

## ACKNOWLEDGEMENTS

We sincerely thank the NASA MEaSUREs Program for providing GRACE-JPL mascon solutions and the Goddard Space Flight Center for providing GSFC mascon solutions. The GRACE-CSR RL05 mascon solutions are provided by the Center for Space Research at the University of Texas at Austin. We thank the Yellow River Conservancy Commission of the Ministry of Water Resources for providing the runoff observations used in this study.

## AUTHOR CONTRIBUTIONS

Y.-P.X. and J.X. designed the study; J.X. did the main calculations and wrote the draft of the manuscript; Y.-P.X. guided the research and revised the manuscript; Y.W. and Y.G. performed data preprocessing.

## CONFLICT OF INTEREST

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## FUNDING

This work was financially supported by the National Key Research and Development Program of China (Grant No. 2018YFC0407401) and the National Natural Science Foundation of China (Grant No. 52009121).

## DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.