Meteorological and Hydrological Drought Risks under Future Climate and Land Use Land Cover Change Scenarios in the Yellow River Basin

1. Introduction

Over recent decades, meteorological and hydrological droughts at regional and global scales have had substantial impacts on agricultural, ecological, and public water systems. Despite rapid scientific and technological advancements, the consequences of such droughts remain unavoidable. For example, approximately 12 billion km² of cropland were affected and over 10 million fatalities occurred due to droughts worldwide between 1900-2010 [1,2]. More recently, severe meteorological and hydrological drought events, including the 2011 East Africa drought, 2012 U.S. Central Great Plains drought, and 2012-2015 California drought, have elicited significant public attention due to their considerable agricultural and societal impacts [3,4,5,6,7]. China has also historically experienced frequent meteorological and hydrological droughts, especially in the ecologically fragile and geographically complex Yellow River Basin (YRB), profoundly influencing agricultural, economic and social development. Numerous studies have examined droughts in the YRB region [8,9,10,11]. For instance, Zhu et al [12] evaluated integrated droughts in the YRB using a comprehensive drought index, finding a drying trend over recent decades, with greater aridity in upstream and northern sections. Li et al. [13] assessed YRB hydrological drought evolution using a nonlinear joint index, identifying increasingly severe hydrological droughts in recent decades, particularly across the Loess Plateau. Wang et al. [14] examined the relationships between hydrological and meteorological droughts in the YRB, again showing worsening hydrological drought conditions. Collectively, these studies highlight the considerable impacts of meteorological and hydrological droughts on the ecologically and socially significant YRB region.

Against the backdrop of global change, meteorological drought patterns in the ecologically significant YRB may exhibit shifts. To understand potential variations under climate change projections, some studies have evaluated future meteorological drought characteristics for the YRB. For instance, Li et al. [15] examined the spatial, temporal, and persistence traits of meteorological droughts under climate warming scenarios across the Loess Plateau region of the YRB. Additionally, Ma et al. [16] assessed meteorological drought risk and change characteristics over a 40-year horizon utilizing eight climate models from the Coupled Model Inter-comparison Project Phase 5 under Representative Concentration Pathway (RCP) 2.6, 4.5 and 8.5 scenarios for the YRB. Research of this nature enhances preparedness for projected meteorological drought shifts in the YRB under global change. Zhang et al. [17] examined the evolution of meteorological drought over the period 2022-2100 under Shared Socioeconomic Pathway (SSP) 2-4.5 and SSP 5-8.5 emissions scenarios across the middle reaches of the YRB, employing copula functions to characterize the joint probability distribution of drought duration and severity. Collectively, these studies exemplify emerging work to anticipate and adapt to potential meteorological drought variations in the YRB resulting from climate change.

However, global change encompasses more than just climate variation; it also includes land use and land cover (LULC) alterations stemming from anthropogenic activities like agricultural expansion, afforestation, and urbanization. LULC characteristics can substantially influence hydrological processes by affecting rainfall infiltration, overland and channel flows, groundwater recharge and discharge, etc. [18]. Modifications to LULC traits may alter regional hydrological cycles and connections between meteorological and hydrological droughts, indirectly exacerbating or relieving hydrological droughts. Therefore, when projecting future hydrological drought risks, it is imperative to account for LULC pattern shifts in addition to climate change. Examining both factors will enable more comprehensive understanding of potential hydrological drought hazards under global change scenarios encompassing climate variations as well as human-induced LULC transformations. This dual perspective is critical for anticipating and managing hydrological drought risks moving forward.

Among the numerous studies on drought risk assessments, the copula function is a popular method that is widely applied in drought risk analysis because of its flexible and margin-free characteristics [19,20,21,22,23]. Zhang et al. [24] employed the copula function to evaluate the meteorological drought risk based on the precipitation data from 1960 to 2005 in the Pearl River basin. Xu et al. [25] used copula function to develop a regional drought risk assessment model, which was applied to analyze the drought frequency in Southwest China during the period of 1961–2012. Chang et al. [26] assessed the multivariate integrated drought risk in the Wei River Basin through the joint probability distribution of drought duration and severity based on Archimedean copulas. Yang et al. [27] assessed the drought risk of the Amu Darya River Basin from 1935 to 2015 by using MCMC-based maximum entropy copula method. However, these studies are mainly focus on historical drought frequency assessments, researches on meteorological and hydrological drought risk assessments under future climate and LULC change scenarios, especially in the YRB, are almost relatively blank.

Therefore, the current study focuses on the ecologically and socially significant YRB, with the primary objectives being: 1) to simulate and analyze projected land use/land cover (LULC) change patterns across the YRB, 2) to construct a drought risk assessment model based on joint probability distribution of drought duration and severity, and 3) to predict and examine meteorological and hydrological drought risks and interconnections under future climate and LULC change scenarios. The remainder of this paper is organized as follows: Section 2 presents the study area and data sources. Section 3 introduces the methodologies utilized in the analysis. Section 4 provides the results and discussion. Finally, Section 5 concludes with a summary of key findings, limitations, and directions for future research. Examining meteorological and hydrological drought projections for the YRB will provide greater understanding of potential changes to this integral system under global change pressures from climate and LULC shifts.

2. Study area and data sources

2.1. The Yellow River Basin

The Yellow River Basin (YRB) is the second largest river basin in China and sixth largest globally. As the cradle of Chinese civilization, it provides water resources for approximately 140 million people while irrigating 15% of China's agricultural land. The YRB spans nine provinces, traversing the Tibetan Plateau, Loess Plateau, and North China Plain before reaching the Bohai Sea. This expansive region exhibits dramatic spatial variability in elevation, precipitation, and temperature. Mean annual precipitation ranges from 123mm to 1021mm, increasing from northwest to southeast. Mean annual temperatures vary from -4°C to 14°C with latitude and elevation [8]. To enable systematic prediction of meteorological and hydrological drought risks across this heterogeneous basin, the YRB was delineated into six subzones based on distinct climatic and topographic characteristics (Figure 1). Detailed traits of each subzone are available in Li et al. [13]. This zonal delineation will facilitate spatial analysis of meteorological and hydrological drought hazards across the climatically and geographically diverse YRB.

2.2. Data sources

2.2.1. NEX-GDDP data

Projections of future climate under varying emissions scenarios are valuable for assessing meteorological and hydrological drought hazards. The National Aeronautics and Space Administration Earth Exchange Global Daily Downscaled Projections (NEX-GDDP) dataset provides statistically downscaled climate projections at high temporal (daily) and spatial (0.25° x 0.25°) resolutions [28,29]. Compared to traditional general circulation models (GCMs), the NEX-GDDP increases spatial resolution and enhances simulation of extreme climate values in the historical period, especially for topographically complex regions. This improves the accuracy of future climate change projections [30]. Several studies utilizing NEX-GDDP for meteorological drought analyses demonstrate its promise as a source of localized and regional-scale climate projections [31,32]. Therefore, this study employs the NEX-GDDP dataset to represent future climate conditions for assessing meteorological drought risk, as well as providing meteorological inputs for evaluating hydrological drought risk. The high-resolution NEX-GDDP data enables robust examination of potential future drought hazards in the topographically and climatically diverse YRB.

Given the BCC-CSM1-1 model within the NEX-GDDP dataset demonstrates skillful simulation and robust predictive capacity for China's climate conditions [33], downscaled projections from this model were selected for the period 2025-2060 under two Representative Concentration Pathway (RCP) scenarios (RCP4.5 and RCP8.5) to evaluate meteorological and hydrological drought risks across the YRB.

2.2.2. Observation data

Meteorological data from 1970-2015, including daily maximum and minimum temperature, precipitation, humidity, solar radiation, and wind speed records from 121 weather stations (Figure 1), were obtained from the National Climate Center of China. To calibrate and validate the SWAT model, naturalized monthly streamflow data from 1970-2015 at five hydrological stations (Figure 1) were acquired from the Yellow River Conservancy Commission. Naturalized streamflow refers to flow that excludes direct anthropogenic influences, such as reservoir storage, irrigation, and industrial, domestic and agricultural water withdrawals [13]. These high-resolution meteorological and hydrological records enable robust calibration and validation of the SWAT model prior to simulating future hydrological conditions across the climatically and topographically diverse YRB.

2.2.3. Spatial data

The spatial data utilized included a digital elevation model (DEM), a digital stream network, a soil map, and four land use/land cover (LULC) maps. The 30m resolution DEM data was extracted from the National Aeronautics and Space Administration's Shuttle Radar Topography Mission (SRTM) Digital Elevation Model. The 1:1,000,000 scale soil map and LULC maps for 1980, 1990, 2000, 2010 and 2020 were acquired from the Data Center for Resources and Environmental Sciences, Chinese Academy of Sciences. This high-resolution spatial data enables accurate characterization of the topographically and environmentally diverse YRB for hydrological modeling and analysis of LULC changes.

3. Methodologies

3.1. Meteorological and hydrological drought indexes

The selection or development of appropriate meteorological and hydrological drought indices is critical for predicting associated risks. Numerous metrics have been proposed for characterizing meteorological and hydrological droughts [13,34,35]. Meteorological drought indices predominantly include the Standardized Precipitation Index (SPI) [36], Precipitation Anomaly Percentage (PAP) [26], Standardized Precipitation Evapotranspiration Index (SPEI) [37], and Reconnaissance Drought Index (RDI) [38]. The SPEI and RDI are well-suited for capturing meteorological drought, as they incorporate precipitation and evapotranspiration, direct drought drivers [39]. Compared to the SPEI, the RDI demonstrates greater sensitivity to environmental changes and clearer evaluation criteria by fitting a distribution function to the precipitation-potential evapotranspiration ratio [40]. Therefore, this study utilizes the RDI as the meteorological drought index. RDI calculation details are available in Tsakiris and Vangelis [41], with classification shown in Table 1.

Several standardized indices have been developed in recent decades for characterizing hydrological drought, including the Surface Water Supply Index (SWSI) [42], Palmer Hydrological Drought Severity Index (PHDI) [43], Standardized Runoff Index (SRI) [44], and Standardized Streamflow Index (SSI) [45]. However, these metrics predominantly characterize surface runoff or streamflow, neglecting the groundwater component. Groundwater (recharge and discharge) also significantly influences hydrological drought as it plays a vital role in hydrological processes [46,47,48]. Therefore, this study employs the Nonlinear Joint Hydrological Drought Index (NJHDI) proposed in our previous study by Li et al. [13] to characterize hydrological drought. The NJHDI comprehensively considers surface water, groundwater, and their nonlinear interconnections. Details on the NJHDI procedures and classification are available in Li et al. [13]. The classification of NJHDI is same as the RDI in Table 1.

3.2. LULC prediction model

3.2.1. CA- Markov model

The Markov chain-cellular automata (CA-Markov) model is a robust approach for modeling and predicting land use/land cover (LULC) change across temporal scales, as it integrates the temporal Markov method with spatial cellular automata techniques [49,50]. The Markov chain is a stochastic modeling method based on the assumption that a system's state probability at a given time can be determined from its earlier state [51]. It has been widely applied to study LULC change dynamics, accurately predicting the magnitude of change at broad scales, but lacks spatial distribution capabilities as spatial categories are not considered [52]. In contrast, the cellular automata (CA) model incorporates spatial attributes by permitting computational experiments on spatial arrangements over time [53]. Coupling CA with Markov modeling has become a predominant technique for LULC change assessment and prediction [54,55,56,57,58].

Therefore, this study utilizes the CA-Markov model to predict future spatial-temporal LULC patterns based on historical trends. The CA-Markov module in IDRISI Andes v17.0 [55] is applied, following five key steps: (1) Calculate LULC transition probability matrix using spatial overlay of historical LULC maps; (2) Determine CA transition rules describing cell state changes; (3) Generate suitability maps defining transition potential; (4) Select 5x5 CA filters to define neighborhood cells; (5) Predict LULC maps for target years based on transition matrix and suitability maps. Further CA-Markov details are available in [54]. This approach couples Markovian temporal modeling with spatial CA techniques to robustly project LULC trends.

3.2.2. Accuracy assessment

In this study, the Kappa statistic is utilized to assess the accuracy and validity of the CA-Markov model. The Kappa statistic serves as a reliable metric for evaluating the agreements in the spatial patterns between two maps [56]. The calculation of Kappa is performed as follows:

$$Kappa=\frac{P_o-P_c}{1-P_c}$$

(1)

where Po represents the proportion of units that agree, while Pc represents the proportion of units in agreement with the expected change. If the value of Kappa falls between 0 and 0.4, it indicates a low level of agreement between the two maps. Conversely, if Kappa ranges from 0.4 to 0.7, it signifies a moderate level of agreement. Finally, if Kappa surpasses 0.7, it suggests a high level of agreement between the two maps.

3.3. SWAT model

The Soil and Water Assessment Tool (SWAT) is a physically based distributed hydrological model that operates in continuous time. It is capable of effectively representing the intricate hydrological processes occurring in watersheds with varying meteorological, land use/land cover (LULC), soil, and geographical conditions, at different temporal scales (daily, monthly, and yearly) [59,60]. Therefore, in this study, the SWAT model is employed to simulate the hydrological processes influenced by climate and LULC changes.

The performance of the SWAT model is assessed using three metrics: the Nash-Sutcliffe coefficient (Ens) [61], coefficient of determination (R²), and percentage bias (PBIAS) [62]. According to Moriasi et al. [63], satisfactory performance of the SWAT model at the monthly scale is indicated by Ens > 0.5, R² > 0.6, and |PBIAS| < 25%.

3.4. Drought event and characteristics definition

The occurrence of a drought event is primarily characterized by drought onset, drought duration, drought severity, and drought termination. Among these, drought duration and drought severity are considered the key indicators for defining a drought event [26]. In this study, the runs theory [64] is utilized to identify a drought event, as well as to determine the duration (D) and severity (S) of the drought. The duration (D) of a drought event is represented by the length of the consecutive dry period, while the severity (S) quantifies the total magnitude of the dry periods, as illustrated in Figure 2.

In Figure 2, two thresholds at -0.5 and 0 are presented on the Y-axis. When the monthly drought index (RDI or NJHDI) falls below -0.5, the corresponding month is classified as a drought month. If the drought index in the subsequent month exceeds 0, or if the drought indexes in the following two months are both within the range of -0.5 and 0, the duration of the drought (D) is denoted as D1, and D1 is equal to the value "a". Alternatively, if the drought index in the second month falls within the range of -0.5 and 0, but then in the third month drops below -0.5, the duration of the drought (D) is referred to as D2, and D2 is calculated as the sum of "b", "d2", and "c". Drought severity (S) is defined as the absolute sum of the drought index values observed during a drought period.

3.5. Copula-based drought risk assessment model

Given that drought duration (D) and drought severity (S) are crucial characteristics of drought events, the assessment of drought risk in this study involves the examination of the joint probability distribution of D and S within a specified range of thresholds. In order to analyze drought risk across various temporal scales, drought events are further categorized into sub-seasonal drought and seasonal drought.

3.5.1. Definition of sub-seasonal drought and seasonal drought

In accordance with the definition of drought event and its characteristics provided in section 3.4, the drought duration (D) represents the count of consecutive dry months within a drought period. On the other hand, the drought severity (S) is quantified as the absolute sum of the drought index values observed throughout the entire duration of the drought, as depicted in equation (2).

$$S=\left|{\displaystyle \sum _{i=1}^n{I}_i}\right|$$

(2)

where I is drought index such as RDI and NJHDI. n is the number of dry months during a drought period.

In the context of a drought event, the duration of the drought is directly proportional to the magnitude of n. To differentiate between short-duration and long-duration drought events, this study introduces two classifications based on the duration of the drought. If the drought duration (D) falls within a single season, i.e., when n is less than or equal to three, the drought event is categorized as a sub-seasonal drought. Conversely, if the drought duration (D) exceeds three months, the drought event is classified as a seasonal drought.

According to the classification of the drought index presented in Table 1, it is established that a month is designated as dry if the corresponding drought index value (I) is less than or equal to -0.5. Furthermore, based on equation (2), it is determined that the severity of the drought (S) during such dry months is greater than or equal to 0.5. In the case of a drought persisting for two or three months, the severity (S) increases to a minimum threshold of 1.0 or 1.5, respectively. Consequently, for a sub-seasonal drought (Figure 3a), the threshold for the severity (S) should be set at least as high as the minimum value among all drought severity (S) in Figure 3a, which is larger than or equal to 0.5. Similarly, for a seasonal drought (Figure 3b), where the drought duration (D) exceeds three months and the corresponding drought severity (S) is greater than or equal to 1.5, the threshold for S should also be set at least as high as the minimum value among all drought severity (S) in Figure 3b, which is larger than or equal to 1.5.

In general, when a drought event is recognized as a sub-seasonal drought or a seasonal drought, the corresponding drought duration (D) and drought severity (S) should be both satisfy the following conditions:

Sub-seasonal drought: $D\le 3,S\ge 0.5$

Seasonal drought: $D\ge 3,S\ge 1.5$

3.5.2. Sub-seasonal and seasonal drought risk assessment models

Copula function can be used to calculate the joint probability distribution of drought duration (D) and drought severity (S) that are both less than or equal to a given threshold as below:

$$P(D\le d,S\le s)=C_{\theta }(F_D(d), F_S(s))$$

(3)

Based on equation (3), the sub-seasonal and seasonal drought risk probability models are established by mathematical deduction as follows:

Sub-seasonal drought risk assessment model:

$$P(D\le d,S\ge s)=F(d)-F(d,s)=F(d)-C_{\theta }(F_{}(d_1), F_{}(s))$$

(4)

Seasonal drought risk assessment model:

$$\begin{array}{l}P(D\ge d,S\ge s)=1-F(d)-F(s)+F(d,s)\\ =1-F(d)-F(s)+C_{\theta }(F_{}(d_{}), F_{}(s))\end{array}$$

(5)

where$F(d)$and$F(s)$are the marginal distribution functions of drought duration (D) and drought severity (S), respectively. $C_{\theta }(F_{}(d_{}), F_{}(s))$ is the copula joint probability distribution.

In this study, we employ the exponential distribution, lognormal distribution, and gamma distribution to fit the marginal distribution of drought duration (D) and drought severity (S). Furthermore, we utilize the Clayton, Frank, and Gumbel-Hougaard copulas to model the joint probability distribution. The selection of the most appropriate marginal distribution and copula functions is determined using the minimum Akaike information criterion (AIC) and bias. These criteria are crucial in ensuring the accuracy and reliability of the selected distributions and copulas.

4. Results and discussion

4.1. LULC prediction in the YRB

The CA-Markov model, which employs equal intervals for simulating and predicting Land Use and Land Cover (LULC) patterns, utilizes the LULC patterns in 1990 and 2000 to predict the LULC in 2010, and the patterns in 2000 and 2010 to predict the LULC in 2020. To assess the accuracy of the CA-Markov model, Kappa coefficients are calculated for 2010 and 2020 by comparing the simulated LULC with observed data using the crosstab module in IDRISI. The resulting Kappa coefficients for 2010 and 2020 are 0.8792 and 0.8850, respectively, both exceeding 0.7, indicating the suitability and reliability of the CA-Markov model for predicting future LULC changes in the YRB.

Based on the CA-Markov model, the historical LULC pattern in 2020 serves as the initial point, combined with historical patterns in 2010, 2000, 1990 and 1980, to predict the LULC patterns in 2030, 2040, 2050 and 2060 using 10, 20, 30 and 40 CA iterations, respectively. The results are illustrated in Figure 4, and quantification of LULC area change characteristics during the period of 2020-2060 is presented in Table 2. Analysis of Figure 4 and Table 2 reveals a significant decline in grassland and cropland areas by approximately 14% in 2060 compared to 2020. Conversely, the areas of the remaining four LULC types are projected to increase. Notably, the build-up area exhibits the most substantial expansion, estimated to be nearly 260% larger in 2060 than in 2020. This expansion can be attributed to rapid urban economic development, increased employment opportunities, and the migration of rural populations to urban areas, thereby accelerating the process of urbanization in the future.

4.2. SWAT model calibration and validation

The YRB is initially divided into 381 sub-basins, ranging in size from 6 km2 to 12,176 km2, based on the Digital Elevation Model (DEM) and channel network. Subsequently, these sub-basins are further partitioned into 1850 hydrological response units (HRUs) using thresholds of 8% for dominant LULC types in 2020, 10% for soil patterns, and 8% for slopes within each sub-basin. This division allows for improved simulation accuracy and a better understanding of the hydrological processes in the YRB. To calibrate and validate the SWAT model for the period of 1970-2015, streamflow observations from the main hydrological stations (see Figure 1) are utilized. These stations, located within each subzone, are selected as they represent the variability and availability of the hydrological processes in the YRB.

This study employs a monthly time step for model runs. The spin-up period utilizes the naturalized stream records from 1970 to generate initial conditions for the model experiment. Model calibration is then performed using the naturalized stream records from 1971 to 1990, while model validation uses data from 1991 to 2015. The calibration and validation results, evaluated using three metrics (i.e., R2 and PBIAS), are presented in Table 3. Analysis of Table 3 reveals that the values of R2 and in both calibration and validation periods exceed 0.6 and 0.7, respectively. Additionally, the absolute values of PBIAS remain within 25%. These results indicate that the SWAT model performs satisfactorily and is capable of reasonably simulating the hydrological processes impacted by climate and LULC changes in the YRB.

4.3. Meteorological drought risk prediction

4.3.1. Selection of the appropriate marginal distribution

In order to determine an appropriate marginal distribution for meteorological drought, the exponential distribution, lognormal distribution, and gamma distribution were employed to model drought duration and severity under RCP 4.5 and RCP 8.5 scenarios. The selection of the best-fitting marginal distributions was based on the minimum values of AIC and Bias, using the AIC and Bias methods. The results are presented in Table 4. It is evident from Table 4 that the error values of the exponential distribution, for both drought duration (D) and drought severity (S), exhibit the lowest values within each subzone under both RCP 4.5 and RCP 8.5 scenarios. This indicates that the goodness-of-fit of the exponential distribution is superior. Therefore, in this study, the exponential distribution is chosen as the appropriate marginal distribution to model meteorological drought duration and severity in the six subzones of the YRB under the two scenarios.

4.3.2. Selection of the appropriate copula function

Once the marginal distribution function is determined, the joint probability of two marginal distributions is modeled using the Clayton, Frank, and Gumbel copulas. The selection of the best-fitting copula is based on the minimum values of AIC and Bias, using the AIC and Bias methods. The fitting errors are presented in Table 5. From the analysis of Table 5, it is evident that the error values of the Frank copula, within each subzone and under both the RCP 4.5 and RCP 8.5 scenarios, are the lowest. This indicates that the Frank copula exhibits the optimal goodness-of-fit effect. Consequently, in this study, the Frank copula is chosen as the appropriate copula function to construct the meteorological drought risk assessment model.

4.3.3. Prediction of meteorological drought risk in the YRB

Utilizing the meteorological drought risk assessment model, predictions are made for the risks of sub-seasonal and seasonal meteorological droughts in the future years of 2025-2060 across the six subzones of the YRB under the RCP 4.5 and RCP 8.5 scenarios. In order to provide a clearer representation of the spatial variation in meteorological drought risk, the spatial risk distributions of sub-seasonal and seasonal meteorological droughts under the two scenarios are illustrated in Figure 5.

Figure 5 reveals that the risks associated with seasonal meteorological droughts remain below 10% throughout the entire YRB under both scenarios. This suggests that the occurrence of seasonal meteorological droughts in the YRB is not expected to be significant in the next 35 years. However, sub-seasonal meteorological drought risks range from 40% to 50% within each subzone, indicating a potential increase in the frequency of sub-seasonal meteorological droughts in the YRB in the future. Additionally, the risk of sub-seasonal meteorological droughts is marginally higher in the upstream subzones (A, B, C) compared to the middle and lower reaches (D, E, F), with risks ranging from 46% to 50% in the upstream and 41% to 45% in the middle and lower reaches. When comparing the RCP 4.5 and RCP 8.5 scenarios, the risk of sub-seasonal meteorological droughts remains relatively stable, except for a 5% increase in the middle stream of zone E under the RCP 8.5 scenario.

4.4. Hydrological drought risk prediction

4.4.1. Selection of the appropriate marginal distribution

In a similar manner to the analysis performed for meteorological drought, the selection of appropriate marginal distributions for hydrological drought is based on the minimum values of AIC and Bias. The results, presented in Table 6, demonstrate the fitting error values for the exponential distribution in relation to both drought duration (D) and drought severity (S) within each subzone under the RCP 4.5 and RCP 8.5 scenarios.

Table 6 clearly indicates that the error values for the exponential distribution, when fitting either drought duration or drought severity, are the lowest within each subzone under both the RCP 4.5 and RCP 8.5 scenarios. This suggests that the exponential distribution provides the best goodness-of-fit. As a result, the exponential distribution is selected as the appropriate marginal distribution in this study to model the hydrological drought duration and severity across the six zones of the YRB under the two scenarios.

4.4.2. Selection of the appropriate copula function

Additionally, in order to select the appropriate copula, the minimum AIC and Bias values are utilized. The fitting errors are presented in Table 7. Analysis of Table 7 reveals that in the PCP 4.5 scenario, the Gumbel copula exhibits the minimum error values in zone A and zone C, while the Frank copula displays the minimum error values in the remaining four zones in the YRB. Consequently, the Gumbel copula is chosen as the suitable copula function for constructing the hydrological drought risk assessment model in zone A and zone C, whereas the Frank copula is selected for constructing the hydrological drought risk assessment model in the other four zones in the YRB. Under the PCP 8.5 scenario, the Frank copula demonstrates the minimum error values in the YRB, with the exception of zone A where the Gumbel copula presents the minimum error value. Therefore, the Gumbel copula is chosen as the appropriate copula function for constructing the hydrological drought risk assessment model in zone A, while the Frank copula is selected for constructing the hydrological drought risk assessment model in the other five zones in the YRB.

4.4.3. Prediction of hydrological drought risk in the YRB

Figure 6 illustrates the spatial distributions of sub-seasonal and seasonal hydrological drought risks under the RCP 4.5 and RCP 8.5 scenarios in the future years of 2025-2060 across the six subzones of the YRB. From Figure 6, it is evident that both sub-seasonal and seasonal hydrological drought risks exhibit significant spatial variability over the next 35 years. Under the RCP 4.5 scenario, sub-zone C and sub-zone E are projected to experience higher risks of sub-seasonal hydrological droughts, ranging from 45% to 50%. However, these two subzones will have lower risks of suffering from seasonal hydrological droughts, as their risks will be less than 10%. Sub-zones A and B will face risks ranging from 20% to 40% for sub-seasonal hydrological droughts, but their risks for seasonal hydrological droughts will decrease by 10% to 20%. Sub-zones D and F will have risks of 10% to 30% for sub-seasonal hydrological droughts, while their risks for seasonal hydrological droughts will increase to 30% to 40%. Under the RCP 8.5 scenario, sub-zone C will be the most susceptible to experiencing sub-seasonal hydrological droughts, with a probability of 45% to 50%. However, it will be the least prone to occur seasonal hydrological droughts, with a probability of less than 10%. Sub-zones A, B, and E will have higher risks of sub-seasonal hydrological droughts, ranging from 20% to 40%. Nevertheless, these subzones will have lower risks of suffering from seasonal hydrological droughts, as their risks will be within the range of 10% to 20%. Conversely, sub-zones D and F will be more susceptible to seasonal hydrological droughts, with probabilities ranging from 30% to 40%. However, their risks for sub-seasonal hydrological droughts will decrease to 20% to 30%.

In general, the YRB is expected to be vulnerable to sub-seasonal hydrological drought, with risks ranging from 30% to 50%. Spatially, zones C and E, which both belong to the Loess Plateau, exhibit the highest risks of experiencing sub-seasonal hydrological droughts, with probabilities ranging from 45% to 50%. It is important to note that the Loess Plateau is characterized by an extremely fragile ecological environment, making recovery from such drought events challenging.

4.5. Analysis of the future concerns regarding sub-seasonal and seasonal droughts

It is challenging to assess whether the sub-seasonal and seasonal droughts, as defined in this study, will be of concern in the future for the YRB without understanding the extent to which future statistics deviate from historical data. Therefore, an analysis of the frequencies of sub-seasonal and seasonal droughts during the historical period from 1970 to 2015 is conducted to gain insights into their potential concerns in the future. The findings of this analysis are presented in Figure 7.

Based on Figure 7, it is evident that the YRB experienced a higher likelihood of sub-seasonal droughts in both meteorology and hydrology compared to seasonal droughts during the historical period. In terms of meteorological drought, the frequency of sub-seasonal drought ranged from 21% to 30% across the YRB, with little variation across different spatial zones. However, the frequency of seasonal drought in the upper reaches of the basin was slightly higher than in the lower reaches. Specifically, zones A and B exhibited drought frequencies of 11% to 20%, while the frequencies in the other four zones were less than 10%. Regarding hydrological drought, there was substantial spatial variation in both sub-seasonal and seasonal drought frequencies during the historical period. Zones B and F had similar frequencies (11% to 20%) for both sub-seasonal and seasonal drought occurrences. Zones C and E had higher probabilities of experiencing sub-seasonal droughts, with frequencies ranging from 21% to 30%. However, these two zones rarely experienced seasonal droughts, with frequencies below 10%. Conversely, zones A and D were more prone to seasonal droughts than sub-seasonal droughts. The frequency of seasonal droughts in these zones was 21% to 30%, while the frequency of sub-seasonal droughts decreased to 0% to 10% in zone A and 11% to 20% in zone D.

Compared to the historical period, the risks of sub-seasonal meteorological and hydrological droughts in the entire YRB will significantly increase in the next 35 years under the RCP 4.5 and RCP 8.5 scenarios. The upstream areas (A, B, C) will experience a 30%-35% higher risk of sub-seasonal meteorological drought under both scenarios, with a range of 46%-50%. In the middle and lower reaches (D, E, F), the risk under the RCP 4.5 scenario will be 25%-30% higher than historical levels, with a range of 41%-45%. However, under the RCP 8.5 scenario, the risks of sub-seasonal meteorological droughts will vary across zone D, zone E, and zone F, with increases of 25%-30%, 30%-35%, and 20% respectively compared to historical levels. Regarding sub-seasonal hydrological droughts, zone C will face the highest risk under both RCP 4.5 and RCP 8.5 scenarios in the future, with a 30%-35% increase compared to historical levels. Similarly, zone E under the RCP 4.5 scenario will experience a 30%-35% higher risk in the next 35 years. However, for the other four zones, the risks of sub-seasonal hydrological droughts will not be a significant concern in the future, as the increases in risk will be below 10% compared to the historical period.

However, in comparison to the historical conditions, the likelihood of experiencing seasonal meteorological and hydrological droughts in the majority of the YRB will remain relatively unchanged in the next 35 years under the RCP 4.5 and RCP 8.5 scenarios. The risks associated with seasonal meteorological drought in the entire YRB under both scenarios are very low, with probabilities below 10%. However, it is important to monitor the downstream area of the YRB (zone F) in the future, as the risk of seasonal hydrological drought in this region will be 20% higher than in the historical period.

In summary, the future sub-seasonal meteorological and hydrological droughts in the YRB should be a cause for concern. Specifically, greater attention should be given to the upstream areas of the YRB (A, B, C) in relation to sub-seasonal meteorological droughts, as the risks in these regions will be 30%-35% higher than in the historical period. Similarly, the Loess Plateau (C, E) should be closely monitored in relation to sub-seasonal hydrological droughts, as the risks in these areas will also be 30%-35% higher than historical conditions. Regarding seasonal droughts, meteorological droughts in the entire YRB will not be a concern in the future, but it is crucial to focus on the downstream area of the YRB (zone F) in the next 35 years, as the risk of seasonal hydrological droughts in this region will increase by 20% compared to historical conditions.

4.6. The relationship between meteorological and hydrological drought risk patterns

Upon analyzing Figure 5 and Figure 6, it becomes evident that significant disparities exist between meteorological and hydrological drought risk patterns. In comparison to hydrological drought, the YRB is projected to be more susceptible to sub-seasonal meteorological drought rather than seasonal meteorological drought in the future, according to two scenarios. This implies that the frequency of meteorological drought will surpass that of hydrological drought. However, the duration of hydrological drought is expected to be longer than that of meteorological drought. The elevated frequency of meteorological drought can be attributed to two primary factors. Firstly, if the intensity of meteorological drought remains relatively mild within a certain time frame, it may not trigger a hydrological drought. Secondly, due to successive meteorological droughts, the short intervals between hydrological droughts make it easier to identify them as drought events. This, in turn, contributes to the higher frequency of meteorological drought. The shorter duration of meteorological drought is primarily a result of the ease with which the balance between precipitation and evaporation can be restored through an increase in precipitation or a decrease in evaporation. However, when it comes to changes in surface water and groundwater, they are more likely to lag behind for extended periods of time. In other words, surface water and groundwater may not recover immediately once the meteorological drought is alleviated. Consequently, the duration of hydrological drought exceeds that of meteorological drought.

From a hydrological cycle perspective, it is widely acknowledged that there exists a definitive linkage between meteorological drought and hydrological drought. The equilibrium between precipitation and evaporation is known to be easily disrupted as meteorological conditions fluctuate, often resulting in the onset of meteorological drought. When meteorological drought persists over an extended period, it engenders two key factors that trigger hydrological drought. Firstly, the prolonged meteorological drought intensifies surface water evaporation and groundwater evaporation, thereby depleting the available surface water and groundwater resources. Secondly, this extended meteorological drought can lead to the thickening and drying of the vadose zone, impeding the runoff yield of surface water and the recharge of groundwater. These two aspects play a pivotal role in inducing hydrological drought, as a prolonged meteorological drought ultimately culminates in the occurrence of hydrological drought. The specific duration of this phenomenon is not only influenced by climate factors but is also closely intertwined with the characteristics of the region's underlying surfaces. Different land use and land cover (LULC) patterns can give rise to distinct mechanisms of runoff generation, consequently exacerbating or mitigating hydrological drought within the same climatic conditions. For instance, woodland and grassland exhibit higher water-holding capacities in the soil, facilitating increased water recharge into groundwater. As a result, when meteorological drought occurs, these LULC patterns can to some extent alleviate the severity of hydrological drought in comparison to urbanized areas.

Based on the land use and land cover (LULC) projections presented in section 4.1, it is evident that urban development, particularly in the middle and downstream regions of the YRB, will experience significant expansion in the future. This expansion of urban areas will have a profound impact on the occurrence of hydrological drought when meteorological drought events transpire. The increased urbanization can hinder the recharge of groundwater due to the presence of impermeable underlying surfaces during periods of precipitation deficit.

This finding further elucidates the reasons behind the projected concerns of sub-seasonal hydrological drought in zone C and zone E, as well as seasonal hydrological drought in zone F, as highlighted in section 4.5. The substantial expansion of urban land in these zones exacerbates the challenges associated with hydrological drought, as the limited recharge of groundwater during meteorological droughts becomes even more difficult due to the impermeable nature of the underlying surfaces.

ident that the error values of the Frank copula, within each subzone and under both the RCP 4.5 and RC

5. Summaries and conclusions

The prediction of both meteorological and hydrological drought risks is essential for establishing early warning systems for drought disasters and formulating effective measures to minimize future drought impacts. In this study, future climate projections for the period 2025-2060 under RCP 4.5 and RCP 8.5 scenarios are obtained using the BCC-CSM1-1 data from the NASA Earth Exchange Global Daily Downscaled Projections (NEX-GDDP) dataset. This downscaled climate dataset provides daily precipitation, as well as daily maximum and minimum temperatures for the YRB. Additionally, the CA-Markov model is employed to simulate and predict future land use and land cover (LULC) change scenarios within the YRB. By combining the future climate scenarios with the projected LULC patterns, a copula-based drought risk assessment model is utilized to estimate the risks of sub-seasonal and seasonal meteorological and hydrological droughts in different regions of the YRB for the next 35 years (2025-2060).

To determine the significance of these future drought events, a comparative analysis is conducted by comparing them to historical drought events. This assessment helps identify whether these projected meteorological and hydrological droughts should be of concern in the future. Lastly, the study analyzes the differences and linkages between the patterns of meteorological and hydrological drought risk, shedding light on their interconnection.

The findings indicate that both sub-seasonal meteorological and hydrological droughts in the YRB will be significant issues in the future. Specifically, greater attention should be given to the sub-seasonal meteorological drought in the upstream areas of the YRB (A, B, C), as the risks in these regions are 30%-35% higher compared to the historical period. Similarly, increased focus is needed on the sub-seasonal hydrological drought in the Loess Plateau (C, E), where the risks are also 30%-35% higher than in the historical period. On the other hand, the future impact of seasonal meteorological drought in the YRB is projected to be negligible, with a low risk level below 10%. However, the downstream area of the YRB (zone F) should be alarmed about the future occurrence of hydrological drought, as the risk will increase by 20% compared to historical conditions. The concerns regarding hydrological droughts in zones C, E, and F primarily stem from the rapid urban land expansion anticipated in the middle and downstream areas of the YRB, which will exacerbate the occurrence of hydrological droughts when meteorological droughts happen.

Great efforts have been made in this study to investigate and predict the spatial risk distributions of sub-seasonal and seasonal meteorological and hydrological droughts scientifically. The results can help decision-makers to find the high-risk occurring zones of meteorological and hydrological droughts more clearly and easily, and hence can prevent and mitigate corresponding local droughts more effectively. In addition, the methods used in this study can apply to other region or basins' for drought risk prediction. Although it has stepped forward an important progress in drought risk prediction with considering not only the variations of future climate conditions, but also the changes in future underlying land surfaces, the present study neglects the influence of future water resource development and utilization (e.g., reservoir storage, irrigation, industrial, domestic and agricultural water consumption) on hydrological drought risk prediction, which will be pursued in our future study.

This study represents a significant effort to scientifically investigate and forecast the spatial distribution of sub-seasonal and seasonal meteorological and hydrological drought risks. The findings have important implications for decision-makers, as they provide clear and accessible information about high-risk areas for meteorological and hydrological droughts. This knowledge can facilitate more effective prevention and mitigation strategies for local drought events. Furthermore, the methodologies employed in this study can be applied to other regions or basins for drought risk prediction.

While this study represents an important advancement in drought risk prediction by considering future climate conditions and changes in underlying land surfaces, it should be noted that the influence of future water resource development and utilization (such as reservoir storage, irrigation, industrial, domestic, and agricultural water consumption) on hydrological drought risk prediction has not been taken into account. This aspect will be addressed in future research endeavors.