Propagation of Meteorological Drought to Agricultural and Hydrological Droughts in the Tropical Lancang–Mekong River Basin

Abstract

In the past several decades, drought events have occurred frequently around the world. However, research on the propagation of drought events has not been adequately explored. This study investigated the drought propagation process from meteorological drought to agricultural drought (PMAD) and from meteorological drought to hydrological drought (PMHD) using a 72-year reanalysis dataset in the tropical Lancang–Mekong River Basin. Firstly, we used a new method—Standardized Drought Analysis Toolbox—to construct drought indices. Then, a linear method (Pearson correlation analysis) and a nonlinear method (mutual information) were used to investigate the drought propagation process. Cross-wavelet analysis and wavelet coherence analysis were employed to explore the statistical relationship among the three drought types. Finally, the random forest method was applied to quantify the major factors in drought response time (DRT). The results revealed the following: (1) both linear and nonlinear methods exhibited strong temporal and spatial consistency for both PMAD and PMHD, with linear relationships being stronger than nonlinear ones. (2) The DRTs of PMAD and PMHD were around 1–2 months and 3–5 months, respectively. Significant differences existed in the DRT between the dry season and the rainy season. (3) A divergent spatial pattern of the proportion of DRT was observed between PMAD and PMHD. (4) Significant statistical correlations between meteorological drought and agricultural drought and between meteorological drought and hydrological drought were observed in specific periods for each sub-region; (5) Hydrometeorological factors contributed the most to DRT, followed by terrain factors and the land cover types. The findings of this study deepened our understanding of the spatial–temporal relationship of multiple drought propagation types in this transboundary river basin.

1. Introduction

Drought is defined as a recurring phenomenon with an abnormal status of hydro-climatic and eco-hydrological elements [1]. As global warming intensifies, meteorological drought events are becoming increasingly frequent, leading to soil moisture drought, hydrological (runoff) drought, and other related issues as a result of the propagation through the water cycle [2,3]. This has led to widespread negative impacts on various aspects, including agriculture [4,5,6,7], water resources [8,9,10], and natural ecosystems [11,12,13], as well as socioeconomic aspects [14,15]. While large negative effects have been ever-lasting, the spatiotemporal patterns of drought propagation and its major influencing factors are largely unknown [16]. Thus, it is of great importance to understand the spatial–temporal patterns of drought propagation and its major influencing factors.

Drought propagation refers to the process of transmission from one drought type to another [17]. Previous research has demonstrated that meteorological droughts are often asynchronous with other types of droughts, meaning that agricultural droughts or hydrological droughts often lag behind meteorological droughts by a certain period [18,19,20]. Propagation between different drought types can be summarized into three categories: (1) The first is the propagation from meteorological drought to agricultural drought (PMAD): the response of soil moisture to precipitation and potential evapotranspiration primarily drives this type of propagation. (2) The second is the propagation from meteorological drought to ecological drought (PMED): the response of vegetation to meteorological drought primarily drives this type of propagation. Vegetation indicators such as Vegetation Condition Index (VCI), Vegetation Health Index (VHI), and Normalized Difference Vegetation Index (NDVI) have been explored and widely used for investigating the impact of meteorological drought on vegetation growth and dynamics [21,22,23,24,25]. (3) The third is the propagation from meteorological drought to hydrological drought (PMHD): the response of hydrological elements (i.e., surface runoff, surface water, groundwater, and streamflow) primarily drives this type of propagation. Previous studies summarized this propagation into four patterns: lagged (presence of a time lag from meteorological drought to hydrological drought), aggregated (multiple meteorological drought events merging into one hydrological drought event), extended (prolonged duration from meteorological drought to hydrological drought), and attenuated (weakening or amelioration of hydrological drought severity due to human activities such as reservoir storage) [16,26]. Though numerous studies have investigated the propagation of various types of droughts mentioned above, the drought response time (DRT) varied with the study region and methodology employed (summarized in

Table 1. The review of previous drought propagation research.

Reference	Study Area	Propagation Type 1	Method 2	Result
[27]	China	PMAD	CA	3.4 months
[28]	Heihe River Basin	PMAD	CA	Average 8 months
[29]	Northeast Asia	PMAD	CA	1–3 months in summer and 5–12 months in winter
[24]	North China Plain	PMED	CA	1.33 months in summer and 2.67 months in winter
[23]	China	PMED	CA	2.67 months in summer and 7 months in winter
[25]	China	PMAD, PMED	CA	1–2.5 months under PMAD and no delay under PMED
[30]	Tarim River Basin	PMHD	CA	2–21 months
[31]	India	PMAD, PMHD	CA	4–5 months under PMAD and 1 month shorter than that under PMHD
[32]	Global	PMAD, PMHD	CA	5.7 months under PMAD and 3.5–14.47 months under PMHD
[33]	Longchuan River Basin	PMAD, PMHD, PHAD	CA	Approximately 2 months under three types of propagation
[34]	Yangtze River Basin	PMAD, PMHD	CA	Less than 2 months under PMAD and 2–6 months under PMHD
[35]	China	PMAD	CA	1–2 months in summer and 2–7 months in the next spring
[36]	Xijiang River Basin	PMHD	RT	Less than 3 months, with the maximum being 78 days
[37]	Yangtze River Basin	PMAD	CA	48 days
[38]	Luanhe River Basin	PMHD	CA	1–7 months in rainy season and 7–12 months in dry season
[39]	Huaihe River Basin	PMHD	RT	1–47 days
[40]	South Korea	PMAD, PMHD	RT	2.83 months under PMAD and 4.34 months under PMHD
[41]	China	PMHD	CA	2–4 months

¹ Propagation Type: PMAD represents the propagation from meteorological drought to agricultural drought, PMHD represents the propagation from meteorological drought to hydrological drought, PMED represents the propagation from meteorological drought to ecological drought, PHAD represents the propagation from hydrological drought to agricultural drought; ² Method: CA represents correlation analysis, RT represents run theory.

).

Generally speaking, the drought propagation process is complicated as it is determined by a lot of factors such as climate conditions, catchment characteristics, and human activities [16]. For example, Zhou et al. [42] investigated the driving mechanisms of PMHD, and their results showed that climate and land use changes were the main factors affecting the propagation process. Han et al. [43] found that soil evaporation and watershed elevation were the main factors influencing the propagation threshold from meteorological drought to groundwater drought in the Xijiang River Basin. Li et al. [29] found that the precipitation and vegetation were major impactors of PMAD in spring, summer, and autumn in Northeast Asia, and the DRT was negatively correlated with precipitation, evapotranspiration, soil moisture and NDVI. Dai et al. [44] investigated the propagation characteristics and mechanism of PMAD and found that an increase in temperature and decreases in soil moisture and precipitation strongly affected the propagation dynamics. In summary, these findings underscore the complexity of understanding and addressing drought propagation, emphasizing the need for comprehensive and integrated strategies in mitigating its impacts.

Extensive studies have been conducted on drought propagation in tropical regions [32,41,45,46,47,48,49]. Yet factors driving the drought propagation in tropical regions were found to be highly complex and divergent for different regions and land cover types. Thus, more research is needed to deepen our understanding of the occurrence and development of drought propagation in tropical regions. As the largest transboundary river basin in Southeast Asia, the Lancang–Mekong River Basin (LMRB) has experienced frequent drought events in recent decades. Quite a lot of studies have assessed the spatiotemporal characteristics of drought in the LMRB [50,51,52,53,54,55]. However, there is a relative scarcity of studies examining the propagation and response of drought, especially for the propagation of multiple drought events. Therefore, there is a need for further in-depth research into the propagation responses of multiple types of drought events in this region. This study aims to investigate the propagation processes of meteorological drought to agricultural drought and hydrological drought in this region and intends to address the following questions: (1) Are there significant differences in the drought propagation time between the two different methods? (2) How long was the DRT of PMAD and PMHD in the study area? Do the two drought propagation types exhibit a significant difference? (3) What factors mainly drive the DRT of the two drought propagation types?

2. Materials and Methods

2.1. Study Area

The tropical Lancang–Mekong River Basin (TRLM) is located in the southeastern part of Asia, encompassing a total area of 795,000 km². Serving as a paradigmatic transboundary area, this region spans various nations’ territories, including Xishuangbanna in China, northeastern Myanmar, eastern Thailand, most regions of Laos and Cambodia, and partial regions of southern Vietnam (Figure 1a). Geographically, the TRLM exhibits a large spatial heterogeneity of altitude pattern with higher elevations in the north, lower elevations in the south, and intermediate elevations in the east and west (Figure 1b). The region is primarily influenced by a tropical monsoon climate with distinct dry and wet seasons. The rainy season lasts from May to October, primarily influenced by the southwestern monsoon, accounting for more than 85% of the total annual rainfall. The dry season persists from November to the next April, predominantly influenced by the northeastern monsoon. During this period, the temperature is lower, and precipitation is relatively scarce. Land cover types in the study region are diverse, with mountainous regions dominated by tropical seasonal rainforests and subtropical evergreen broad-leaved forests and plain regions consisting of farmland, orchards, wetlands, shrubs, bodies of water, and artificial grasslands (Figure 1c). Based on the spatial heterogeneity of topography and vegetation distribution, we further divided the study region into six sub-regions from north to south: Lower Lancang River, Northern Highlands, Korat Plateau, Annan Mountains, Tonle Sap Lake Basin, and Mekong Delta Region (Figure 1a).

2.2. Data Description

2.2.1. Meteorological Reanalysis Data

The meteorological data were sourced from the ERA5-Land reanalysis dataset (

Table 2. Summary of the dataset used in this study.

Data Type	Spatiotemporal Resolution	Time Span	Data Source
ERA5-Land reanalysis data	0.1° (~10 km), monthly	1950–2021	https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-land-monthly-means?tab=overview (accessed on 4 October 2022)
GlobeLand30 land cover	30 m, yearly	2020	httes://www.globallandcover.com (accessed on 12 July 2023)
GIMMS_3g NDVI	1/12° (~8 km), 15 days	1981–2015	https://ecocast.arc.nasa.gov/data/pub/gimms/ (accessed on 27 July 2023)
SRTM DEM	30 m	2000	https://earthexplorer.usgs.gov/ (accessed on 8 September 2022)

). Compared to the ERA5 dataset, ERA5-Land demonstrates a superior capability in representing elements of the water cycle [56]. In addition, previous research had indicated that ERA5-Land exhibited a higher accuracy in comparison to most other reanalysis datasets (i.e., MERRA-2, JRA-55, GLDAS) [57,58,59]. In this study, precipitation, soil moisture, and runoff data from the ERA5-Land dataset were employed to calculate three drought indices, namely the Standardized Precipitation Index (SPI), Standardized Soil Moisture Index (SSMI), and Standardized Runoff Index (SRI). The data had a monthly temporal resolution and a spatial resolution of 0.1° (~10 km), covering the period from 1950 to 2021.

2.2.2. Land Cover Data

The land cover data originated from the GlobeLand30 dataset (Table 2). This dataset is a globally comprehensive land cover dataset with a spatial resolution of 30 m [60]. GlobeLand30 was mapped using a pixel–object–knowledge-based classification algorithm based on multiple satellites, including Landsat TM5, ETM+, OLI, and Chinese HJ-1 satellite images. The 2020 version used in this study also incorporated high-resolution multispectral images from the GaoFen-1 satellite with a resolution of 16 m. With an overall accuracy of more than 85% globally, this dataset has been extensively utilized in global and basin-scale studies [61,62,63,64,65].

2.2.3. Vegetation and Terrain Data

The vegetation dataset was derived from the GIMMS NDVI3g dataset (Table 2), which features a temporal resolution of 15 days and a spatial resolution of 1/12° (~8 km). Covering the period from 1981 to 2015, it is one of the longest-existing NDVI remote sensing datasets [66]. This dataset has been extensively employed in the assessment of long-term trends in vegetation activity [67,68,69]. The terrain dataset utilized in this study was obtained from the Shuttle Radar Topography Mission (SRTM) (accessed by https://earthexplorer.usgs.gov/ (accessed on 8 September 2022)). The spatial resolution of the SRTM is 30 m. The vegetation and terrain data were employed in conducting the driving factor analysis of DRT.

2.3. Methods

2.3.1. Standardized Drought Indices

The SPI, SSMI, and SRI were employed to characterize meteorological drought, agricultural drought, and hydrological drought, respectively. The conventional methodology of developing these three drought indices often relies on comparing representative parameter probability distribution functions. However, such probability distribution functions often lack spatial universality, making it difficult to apply them in continental- or global-scale studies [70]. The Standardized Drought Analysis Toolbox (SDAT), developed by Farahmand and AghaKouchak [71], is a new method for deriving non-parametric univariate and multivariate standardized drought indices. This approach does not require parameter estimation and goodness-of-fit evaluation, which makes it much more efficient when conducting large-scale (continental or global) studies [71]. In this study, we used SDAT to construct time series of SPI, SSMI, and SRI drought indices spanning from 1950 to 2021. The formulations involved are outlined as follows:

$$p\left(x_i\right)=\frac{i-0.44}{n+0.12}$$

(1)

$$SI=\left\{\begin{array}{c}-\left(t-\frac{c_0+c_{1}t+c_2{t}^2}{1+d_{1}t+d{t}^2+d_3{t}^3}\right), 0<p\le 0.5\\ +\left(t-\frac{c_0+c_{1}t+c_2{t}^2}{1+d_{1}t+d{t}^2+d_3{t}^3}\right), 0.5<p\le 1\end{array}\right.$$

(2)

$$t=\left\{\begin{array}{c}\sqrt{ln\frac{1}{p^2}}, 0<p\le 0.5\\ \sqrt{ln\frac{1}{{\left(1-p\right)}^2}}, 0.5<p\le 1\end{array}\right.$$

(3)

where SI is the standardized drought index, x is the time series data for hydrometeorological variables, n is the sample size, i denotes the rank of non-zero precipitation (or soil moisture or runoff) data from the minimum, and $p\left(x_i\right)$ is the corresponding empirical probability. As for the constants, $c_0$ = 2.515517, $c_1$ = 0.802583, $c_2$ = 0.010328, $d_1$ = 1.432788, $d_2$ = 0.189269, and $d_3$ = 0.001308.

In this study, SPI, SSMI, and SRI were calculated at a 1-month time scale using the Equations (1)–(3). Namely, we also defined SPI-n using the first n months’ hydrometeorological variables to calculate SPI. For example, SPI-3 denotes a 3-month precipitation accumulation period.

2.3.2. Determination of DRT

Previous studies have demonstrated that both linear and nonlinear relationships jointly govern the drought propagation process [72]. Therefore, we simultaneously utilized linear and nonlinear methods to investigate the PMAD and PMHD. The Pearson correlation coefficient (PCC) has been widely employed as a fundamental tool for elucidating linear relationships between different drought types and can effectively show the response intensity and time [73]. The integration of mutual information, originating from information theory and entropy theory, was adopted in capturing nonlinear relationships [74,75].

• The linear method—Pearson correlation coefficient

The cetermined by the time scales with the largest PCC calculated for each SPI-n (n = 1, 2,..., 12) and SSMI. Similarly, the linear DRT of PMHD was determined by the time scales with the largest PCC calculated for each SPI-n (n = 1, 2,..., 12) and SRI. The specific calculation process is as follows:

$$PC{C}_{X,Y}=\frac{cov\left(X,Y\right)}{{\sigma }_X{\sigma }_Y}=\frac{E\left(XY\right)-E\left(X\right)E\left(Y\right)}{\sqrt{E(X^2)-E^2\left(X\right)}\sqrt{E(Y^2)-E^2\left(Y\right)}}$$

(4)

$$r_n=PC{C}_{SP{I}_n,SSMI or SRI}, 1\le n\le 12$$

(5)

$$MPCC=\max \left(r_n\right), 1\le n\le 12$$

(6)

$$DR{T}_{linear}=N \left(r_N=MPCC\right)$$

(7)

where $r_n$ represents the PCC between SPI and SSMI or SRI for a given time scale n ($1\le n\le 12$), MPCC is the maximum value among the 12 $r_n$ values, and the corresponding time scale n associated with the MPCC indicates the linear DRT of PMAD or PMHD. The PCC value ranges between −1 and 1, where a larger absolute value indicates a stronger correlation.

2.

The nonlinear method—mutual information

In the process of drought propagation and evolution, various meteorological factors and basin characteristics can collectively influence the spread of drought, potentially leading to nonlinear relationships that cannot be adequately addressed through Pearson correlation analysis. To address this limitation, mutual information (MI) was used to quantify the nonlinear correlations between meteorological drought and agricultural or hydrological drought. The calculation of MI is as follows:

$$M{I}_{X,Y}={{\displaystyle \sum }}_Y{{\displaystyle \sum }}_{X}p\left(x,y\right)\log \left(\frac{p\left(x,y\right)}{p\left(x\right)p\left(y\right)}\right)$$

(8)

$$r_n=M{I}_{SP{I}_n,SSMI or SRI}, 1\le n\le 12$$

(9)

$$MMI=\max \left(r_n\right), 1\le n\le 12$$

(10)

$$DR{T}_{non-linear}=N \left(r_N=MMI\right)$$

(11)

where $p\left(x,y\right)$ is the joint probability distribution function of X and Y; $p\left(x\right)$ and $p\left(y\right)$ are the marginal probability density functions of X and Y respectively; $r_n$ represents the MI between SPI and SSMI for PMAD or between SPI and SRI for PMHD at a given time scale n ($1\le n\le 12$); MMI is the maximum value among the 12 $r_n$ values; and the corresponding time scale n associated with the MMI indicates the nonlinear DRT of PMAD or PMHD.

2.3.3. Determination of the Relationship between Different Types of Droughts

This study employed cross-wavelet transform (XWT) and wavelet coherence (WTC) to explore the relationships between meteorological drought and agricultural drought or hydrological drought. Cross-wavelet analysis is one of the most effective tools based on cross-spectral analysis and wavelet analysis to investigate the correlation between two time series [76,77,78]. It builds upon continuous wavelet analysis to analyze significant periods, corresponding time spans, phase relationships, and time lags between two sets of time series. XWT can elucidate coherence and phase relationships between different time series, while it may lack precision in terms of consistency in some conditions. WTC overcomes this limitation, providing accurate results about coherence relationships [79]. The range of WTC values is from 0 to 1, where values closer to 1 indicate stronger correlations between two time series. Both of these methods have been extensively employed in various meteorological and hydrological studies [19,79,80]. In this study, the MATLAB r2022b wavelet analysis toolbox (accessed by https://atoc.colorado.edu/research/wavelets on 15 October 2022) developed by Grinsted et al. [76] was used.

2.3.4. Variable Importance Based on Random Forest

The random forest method has been widely used to measure the variable importance in previous research [81,82]. It was also employed in this study to quantify the importance of nine hydrometeorological and environmental factors: precipitation, runoff, soil moisture, temperature, elevation, slope, and proportions of three major land cover types (i.e., forest, shrub, and crop). The sum of all the importance values is 1, and a larger importance value indicates that a variable plays a more important role in affecting drought propagation.

3. Results

3.1. Spatial and Temporal Characteristics of the Three Drought Types

Figure 2 shows the Theil–Sen linear regression of each drought index. The SPI, SSMI, and SRI exhibited a high degree of spatial and temporal consistency in their trends, showing a distinct spatial heterogeneity with an insignificant decreasing trend dominating in most parts of the upper region and an insignificant upward trend dominating in most parts of the lower region. Specifically, regions with insignificant upward trend of three drought indices included the southwestern part of the Korat Plateau, the southern part of the Annan Mountains, the southeastern and southwestern parts of the Tonle Sap Lake Basin, and most areas of the Mekong Delta. This suggested a lower likelihood of these regions experiencing drought events. Regions with a decreasing trend of three drought indices included the Lower Lancang River, the Northern Highlands, and most portions of the Korat Plateau. Notably, the downward trend was more pronounced in the Lower Lancang River and the northeastern part of the Korat Plateau, implying a drier climate and a higher probability of drought events occurring in these regions. Across most of the study area, trends of the three drought indices were relatively insignificant.

3.2. The Process of PMAD and PMHD

The results from both linear correlation and nonlinear methods indicated a significant spatial heterogeneity of the correlation relationship between meteorological and agricultural droughts and between meteorological and hydrological droughts (Figure 3). Taking the PMAD for example, the spatial pattern revealed high MPCC and MMI values in plateaus, lake basins, and deltas, while mountainous regions exhibited relatively low values. For instance, the MPCC and MMI values between meteorological and agricultural droughts exceeded 0.75 and 0.50, respectively, in sub-regions such as the Korat Plateau, Tonle Sap Lake Basin, and Mekong Delta. Conversely, the MPCC and MMI values in regions like the Lower Lancang River, the Northern Highlands, and the southern part of the Annan Mountains were below 0.70 and 0.40, respectively (Figure 4(a1,a2)). As for PMHD, similarly, the spatial pattern revealed high linear correlation and MI values in lake basins and deltas, while plateau regions exhibited relatively low values. Furthermore, it can be observed that results from both linear and nonlinear methods had good spatial consistency for both PMAD and PMHD. Regions with higher MPCC values also tended to have higher MMI values (Figure 3(a3,c3)).

The DRT of PMAD and PMHD produced by both linear and nonlinear methods also exhibited strong spatial consistency (Figure 3(b1–b3,d1–d3)). More than 80% of the study area showed identical DRTs of PMAD and more than 70% of the study area showed identical DRTs of PMHD under the two methods. As for PMAD, the remaining 20% of the region exhibited less DRT from the linear method than from the nonlinear method, indicating a potentially more complex agricultural drought propagation process in these areas (Figure 3(b3)). For example, a significant difference between the linear and nonlinear propagation times can be observed in the Annan Mountains and the Northern Highlands. In the Annan Mountains, most areas showed 1 month of DRT from the linear method (Figure 3(b1)), while the area proportions of 1 month and 2 months were approximately the same for the DRT from the nonlinear method (Figure 3(b2)). Considering the PMHD, the remaining 30% of the region exhibited various DRT differences, with 14% of the region exhibiting less DRT from the linear than nonlinear method, and 16% of the region showed the opposite distribution, indicating a potentially more complex hydrological drought propagation process in these areas (Figure 3(d3)). For example, a significant difference between the linear and nonlinear DRTs can be observed in the Northern Highlands and a large part of the Korat Plateau and Annan Mountains.

Both methods revealed that the DRT was between 1 and 2 months of PMAD in the study area and 1 and 8 months of PMHD in most areas except for the Korat Plateau, for which the DRT ranges between 1 and 12 months (Figure 3(b1,b2,d1,d2)). Areas with an agricultural DRT of 2 months were mainly located in the Lower Lancang River, Northern Highlands, and Korat Plateau, while areas with a DRT of 1 month were predominantly situated in the southern part of the Annan Mountains, Tonle Sap Lake Basin, and Mekong Delta (Figure 3(b1,b2)). Meanwhile, areas with a hydrological DRT of 1–3 months were mainly located in the Lower Lancang River, Northern Highlands, and central Annan Mountains, while areas with a DRT of 4–6 months were predominantly situated in most of the Tonle Sap Lake Basin and Mekong Delta region. It is worth noting that the Korat Plateau exhibited remarkable spatiotemporal differentiation with a hydrological DRT of 1–12 months, indicating an extremely complex hydrological drought propagation process in this region (Figure 3(d1,d2)).

The proportion of the DRT in each sub-region shown in

Table 3. The DRT proportion in each sub-region under PMAD.

Sub-Region	Method	1 Month	2 Months	3 Months or More
Sub-region I	linear	33.02%	66.98%	0
	nonlinear	15.41%	84.59%	0
Sub-region Ⅱ	linear	39.91%	60.09%	0
	nonlinear	15.41%	84.59%	0
Sub-region Ⅲ	linear	41.39%	58.61%	0
	nonlinear	32.96%	67.04%	0
Sub-region Ⅳ	linear	83.33%	16.67%	0
	nonlinear	50.62%	49.38%	0
Sub-region Ⅴ	linear	64.07%	34.91%	1.02%
	nonlinear	49.03%	49.77%	1.2%
Sub-region Ⅵ	linear	71.39%	28.61%	0
	nonlinear	66.67%	33.33%	0

and

Table 4. The DRT proportion in each sub-region under PMHD.

Sub-Region	Method	1–3 Months	4–6 Months	7–9 Months	10–12 Months
Sub-region Ⅰ	linear	70.76%	29.24%	0	0
	nonlinear	76.42%	23.58%	0	0
Sub-region Ⅱ	linear	80.48%	18.91%	0.61%	0
	nonlinear	78.26%	21.07%	0.67%	0
Sub-region Ⅲ	linear	15.26%	50.11%	28.91%	5.72%
	nonlinear	12.18%	49.73%	28.83%	9.26%
Sub-region Ⅳ	linear	68.91%	30.25%	0.84%	0
	nonlinear	62.92%	36.67%	0.41%	0
Sub-region Ⅴ	linear	36.13%	53.25%	9.61%	1.01
	nonlinear	35.34%	51.26%	9.91%	3.49%
Sub-region Ⅵ	linear	32.77%	64.84%	1.68%	0.71%
	nonlinear	33.26%	62.21%	1.92%	2.61%

indicated that notable differences were found between the proportion of the agricultural DRT while minor differences existed in the proportion of hydrological DRT for the linear method and that for the nonlinear method. Taking the DRT of PMAD for example (Table 3), the Northern Highlands and the Annan Mountains exhibited the greatest proportion gaps between the linear and nonlinear DRTs, with 24.5% and 32.71%, respectively. The smallest difference in the DRT proportion between the two methods was observed in the Korat Plateau and the Mekong Delta, with 8.43% and 4.72%, respectively. Significant differences can also be observed between the proportion of 1 month and 2 months in most regions. In general, the agricultural DRT proportion of 2 months was larger than that of 1 month in the upper sub-regions (i.e., the Lower Lancang River, the Northern Highlands, and the Korat Plateau) while the proportion of 2 months was smaller than that of 1 month in the lower sub-regions (i.e., the Annan Mountains, Tonle Sap Lake Basin, and Mekong Delta). As for the DRT of PMHD (Table 4), the largest proportion gap between the linear and nonlinear DRTs was 6.4% in the Annan Mountains with a DRT of 4–6 months. Significant differences can also be observed between the proportion of various DRTs in most regions. For example, the proportion of 1–3 months was larger than that of 4–6 months in the Lower Lancang River, the Northern Highlands, and the Annan Mountains, while it was smaller than that of 4–6 months in the Korat Plateau, Tonle Sap Lake Basin, and Mekong Delta.

Figure 5 displays the DRTs of PMAD and PMHD in the dry and rainy seasons for two different methods. The DRT in most sub-regions was relatively consistent under the two methods. Taking the DRT of PMAD into consideration (Figure 5(a1,a2)), it is evident that the DRT in the rainy season was generally 1 month across the six sub-regions. However, during the dry season, the DRT was generally longer, predominantly 2 months. For most sub-regions, the linear and nonlinear DRTs were consistent for two seasons. For instance, in the Lower Lancang River, Northern Highlands, Korat Plateau, and the Annan Mountains, the linear and nonlinear DRTs were 1 month during the rainy season and 2 months during the dry season. There are also notable differences between the two methods in some regions. For example, in the Tonle Sap Lake Basin and the Mekong Delta, the linear DRTs were both 1 month during the dry and rainy seasons, while the nonlinear DRTs were 1 month during the dry season and 2 months during the rainy season. This suggests that the process of PMAD was more complex in these two sub-regions during both dry and rainy seasons. Additionally, both the MPCC and MMI values of the dry season were substantially greater than those of the rainy season. This indicates that meteorological and agricultural droughts were more strongly linked during the dry season. As for the DRT of PMHD (Figure 5(b1,b2)), it is evident that the DRT of the dry season was generally 1 month across the six sub-regions, while it was generally longer during the rainy season, predominantly 3–5 months. For most sub-regions, the linear and nonlinear DRTs were consistent for both seasons. For instance, in the Lower Lancang River, Northern Highlands, Annan Mountains, and the Tonle Sap Lake Basin, the linear and nonlinear DRTs were both 4 months (except for Northern Highlands with a 3-month DRT) during the rainy season and 1 month during the dry season. Similar to the PMAD, both the MPCC and MI values during the dry season were notably greater than those during the rainy season. This also indicates that meteorological and hydrological droughts were more strongly linked during the dry season.

Given that both the linear and nonlinear methods showed the effective estimation of DRT in the study area, the final DRT of each sub-region was determined by the maximum MPCC and MI values among all 12 MPCC and MI values separately calculated for SPI-n and the one-month scaled SSMI (for PMAD) or SRI (for PMHD) (

Table 5. The DRT in each sub-region under PMAD and PMHD.

Sub-Region	PMAD		PMHD
	Linear Method	Nonlinear Method	Linear Method	Nonlinear Method
Sub-region Ⅰ	2 months	2 months	3 months	3 months
Sub-region Ⅱ	2 months	2 months	3 months	3 months
Sub-region Ⅲ	2 months	2 months	5 months	5 months
Sub-region Ⅳ	1 month	1 month	3 months	3 months
Sub-region Ⅴ	1 month	2 months	4 months	4 months
Sub-region Ⅵ	1 month	1 month	4 months	4 months

). The results showed that the linear DRTs of PMAD were 2 months in the Lower Lancang River, Northern Highlands, and Korat Plateau and 1 month in the Annan Mountains, Tonle Sap Lake Basin, and Mekong Delta Region. Correspondingly, the nonlinear DRTs in the above-mentioned sub-regions were 2 months and 1 month, respectively. It is obvious that the linear DRT and nonlinear DRT are the same almost in all the sub-regions except for the Tonle Sap Lake Basin, where the former was shorter than the latter. Given that the MPCC was much larger than the MMI value in the Tonle Sap Lake Basin, we suggest that the linear method is more suitable for predicting the process of the PMAD. It is worth noting that the linear and nonlinear DRTs of PMHD were consistent in all sub-regions. Therefore, the final DRTs of PMHD were 3 months, 3 months, 5 months, 3 months, 4 months, and 4 months in the six sub-regions sequentially.

3.3. The Relationship between Meteorological Drought and Other Two Drought Types

The cross-wavelet power spectra (XWT) and wavelet coherence (WTC) analyses were conducted based on the SPI-n (n was determined based on DRT) and the SSMI or SRI of each sub-region. These two methods allow us to investigate the relationships between meteorological drought and agricultural drought, as well as between meteorological drought and hydrological drought.

3.3.1. The Relationship between Meteorological Drought and Agricultural Drought

The results of XWT and WTC showed that significant statistical correlations between agricultural drought and meteorological drought can be observed for each sub-region, but the significant periods for each sub-region were different (Figure 6(a1–f1)). Specifically, for the Lower Lancang River, significant correlations can be observed mainly during 1959–1970 and 1972–1980, with periods of 28–40 months and 48–72 months, respectively. The Northern Highlands showed significant correlations primarily during 1955–1966, with a period of 40–60 months. The Korat Plateau sub-region exhibited significant correlations during 1965–1968, 1975–1980, and 1995–2005, with periods approximately ranging from 8 to 48 months. The Tonle Sap Lake Basin exhibited significant correlations during 1970–1984 and 1995–2005, with a period of 16–48 months. The Annan Mountains exhibited significant correlations during 1995–2005, with a period of 16–64 months. The Mekong Delta exhibited significant correlations during 1992–2006, with a period of 12–64 months. It is noteworthy that all the sub-regions exhibited a substantial amount of short-period wavelet power. Furthermore, the cross-wavelet power duration of the low-frequency range between SPI-n and SSMI was relatively short at the 95% confidence level, indicating an unstable relationship between meteorological and agricultural drought. It should be noted that both SPI-n and SSMI showed significant co-variability at each time scale. The WTC plots (Figure 6(a2–f2)) indicated that SPI-n and SSMI exhibited significant coherency in most parts of each sub-region. This demonstrated that meteorological drought was the main driver of agricultural drought in this study region.

3.3.2. The Relationship between Meteorological Drought and Hydrological Drought

The results of XWT and WTC showed that significant statistical correlations between hydrological drought and meteorological drought can be observed for each sub-region, but the significant periods for each sub-region were different (Figure 7(a1–f1)). For example, in the Lower Lancang River, significant correlations can be observed mainly during 1952–1956 and 1966–1972, with periods of 4–16 months. The Northern Highlands showed significant correlations primarily during 1952–1956, with a period of 6–16 months. The Korat Plateau sub-region exhibited significant correlations during 1965–1968 and 1998–2002, with periods approximately ranging from 12 to 32 and 32 to 48 months, respectively. The Tonle Sap Lake Basin exhibited significant correlations during 1974–1982 and 1995–2002, with a period of 32–48 months. The Annan Mountains exhibited significant correlations during 1960–2010 and 1995–2005, with periods of 64–160 and 32–64 months, respectively. The Mekong Delta exhibited significant correlations during 1992–2006 and 1980–2010, with periods of 20–48 and 112–160 months, respectively. Note that all the sub-regions exhibited a substantial amount of short-period wavelet power. Furthermore, the cross-wavelet power duration of the low-frequency range between SPI-n and SRI was relatively short at the 95% confidence level, indicating an unstable relationship between meteorological and hydrological drought. Both SPI-n and SRI also showed significant co-variability at each time scale. The WTC plots (Figure 7(a2–f2)) indicated that SPI-n and SRI exhibited significant coherency in most parts of each sub-region, implying that meteorological drought was a key driver of hydrological drought.

4. Discussion

4.1. Reasons for the DRT Variations in Each Sub-Region

Figure 8 presents the influence of hydrometeorological and environmental factors on the DRT of PMAD. It was found that higher precipitation, soil moisture, and runoff typically corresponded to shorter agricultural DRTs. This may imply that wetter regions exhibit shorter DRTs of PMAD. In contrast to the above hydrometeorological factors, areas with slightly higher elevation presented longer agricultural DRTs than those with lower elevation. This may be attributed to the high-elevation areas which generally have a faster surface runoff and soil water movement. We also found temperature and slope exhibited a small difference in different DRTs. This implies that these two factors had subtle influences on the DRT of PMAD. Meanwhile, the land cover type also served as an important factor, with a low proportion of cropland and a high proportion of forest tending to have shorter DRTs of PMAD. This may be attributed to forest-dominated regions often having higher evapotranspiration and the water resources being consumed more, thus leading to a relatively shorter DRT [83], while in cropland-dominated regions irrigation measures are often adopted to stabilize the soil moisture when faced with meteorological drought events, further prolonging the propagation process.

Figure 9 presents the influence of hydrometeorological factors and basin characteristics on the DRT of PMHD. Aligning with the findings for PMAD, regions with abundant precipitation, runoff, and soil moisture tended to have shorter hydrological DRTs. This also implies that wetter regions exhibited shorter DRTs of PMHD. Compared with water-related factors, temperature presented a contrasting result: regions with higher temperatures showed slightly longer DRTs of PMHD. Regions with higher elevations and slopes tended to exhibit shorter DRTs of PMHD, implying that topography played a crucial role in affecting the DRT of PMHD. Similar to PMAD, land cover types displayed a consistent result: areas with a low proportion of cropland and a high proportion of forest tend to have a shorter DRT of PMHD, indicating that land cover type significantly influences the hydrological propagation process.

The hydrometeorological factors were further quantified by using the random forest method. As illustrated in Figure 10, apart from the factors that directly influence the propagation process (soil moisture for PMAD and runoff for PMHD), precipitation was the component that contributed the most to DRT. This is reasonable as precipitation is the major source of soil moisture and surface runoff [84]. Overall, the meteorological factors contributed the most to DRT, followed by terrain factors, and the land cover types contributed the least. It is worth noting that the terrain factors were more important than land cover types for both drought propagation types.

4.2. Comparison with Previous Research

In this study, it is observed that the DRT under PMAD and PMHD was 1–2 months and 1–9 months, respectively, in most areas (Figure 3(b1,b2,d1,d2)). The DRT exhibited a spatiotemporal distribution pattern of shorter in the north and longer in the south. Meanwhile, among the six sub-regions, sub-region Ⅲ (Korat Plateau) displayed significant spatiotemporal differentiation, which is consistent with previous research [55,85,86]. The results also indicated that the DRT was shorter during the dry season and longer during the rainy season, which was highly consistent with previous studies [41,45,84]. However, inconsistent results can also be found in some other previous studies. For example, Wang et al. [87] explored the propagation process of meteorological and hydrological drought to socioeconomic drought. The result showed the DRT and correlation of meteorological and hydrological drought with socioeconomic drought in the wet season were higher than those in the dry season. Xu et al. [38] found the DRT of PMHD was shorter in the rainy season (1–7 months) and longer in the dry season (7–12 months) in the Luan River Basin. Zhu et al. [88] also showed a similar pattern. Overall, the factors affecting drought propagation are complex.

4.3. Uncertainties, Limitations, and Future Direction

The hydrometeorological elements from the ERA5-Land reanalysis dataset inevitably introduced certain errors compared to ground-based measurements. However, ERA5-Land offers a higher spatiotemporal resolution and better hydrological representation compared with other remote sensing and reanalysis datasets, making it widely applied in global and regional drought research [89,90,91,92]. For instance, Dergunov et al. [93] compared the temperature data of ERA5-Land, Climate Research Unit gridded time series (CRU TS), and the Global Forecast System (GFS) and found that although the three products had a high correlation with temperature, the ERA5-Land dataset was more suitable for capturing small changes due to its higher spatiotemporal resolution. The evaluation of soil moisture and surface runoff data by previous studies indicated the outperformance of the ERA5 product [94]. Specifically, Li et al. [94] evaluated the soil moisture from 25 networks of situ observations and 5 reanalysis products (IRA-55, CFSR, ERA-Interim, ERA5, and MERRA-2). The result showed that the ERA5 soil moisture had a higher consistency and closer intensity. Harrigan et al. [95] evaluated the surface runoff data of ERA5-Land and the Global Flood Awareness System (GloFAS) on a global scale. The findings indicated that runoff from ERA5-land was highly credible in 86% of the global discharge stations. Additionally, Bain et al. [96] assessed three surface runoff products from ECMWF (i.e., ERA-Interim, ERA5, and ERA5-Land) in the Mississippi River Basin and found that ERA5-Land had the smallest runoff data bias and the highest correlation with model-simulated runoff. In terms of precipitation, we compared the applicability of the ERA5-Land and CRU TS 4.06 datasets with gauged precipitation data in the study area (Figure A1). The results showed that the root mean square error (RMSE) and R² of ERA5-Land precipitation were 23% lower and 0.24 higher than those from CRU TS, respectively, indicating the superiority of ERA5-Land over the CRU TS product.

Current widely used methods for studying drought propagation mainly include correlation analysis, run theory, and nonlinear response function methods. The advantage of run theory lies in its applicability to different drought indices based on a drought threshold. However, agricultural and hydrological drought events may respond to meteorological drought at different time scales (i.e., monthly, seasonal, yearly), which would potentially lead to significant discrepancies when based on fixed-time-scale drought indices. Additionally, the determination of drought thresholds lacks a unified protocol. This often leads to different drought propagation results for the same research [16,97]. The nonlinear response function method can simulate the nonlinearity in drought propagation processes. Thus, it can reveal the internal connections among different types of droughts. However, the construction of response functions is subjective and lacks generality. Correlation analysis simultaneously accounts for both linear and nonlinear relationships between two objective drought types. It also considers the entire time series cycle of wet and dry conditions. Moreover, it is simple in calculation and is widely applicable (Table 1).

This study investigates the propagation characteristics of drought based on monthly meteorological, agricultural, and hydrological drought indices. However, it is important to acknowledge that drought propagation may occur at finer temporal scales. The occurrence and spread of agricultural drought, for instance, may take place on the weekly scale. In such cases, an analysis based on monthly agricultural drought indices may not adequately capture the response of crops during specific critical periods. Therefore, future research can employ drought indices with a finer temporal resolution [16]. Moreover, due to the lack of 72-year land cover data, we only used the land cover data for 2020. This would undoubtedly bring some uncertainty to this study, but the uncertainty is limited. The reason is that the proportion of land cover change in this basin was small (5.27% between 1997 and 2010) in recent decades according to Sam et al. [98]. Future research needs to concentrate more on the impact of land cover change on drought propagation. Additionally, this study did not explore the process of propagation occurrence and development. Subsequent research could enhance the understanding of the evolution of drought by employing relevant hydrological models.

5. Conclusions

This study investigated the propagation process from meteorological drought to agricultural drought (PMAD) and from meteorological drought to hydrological drought (PMHD) in the TRLM based on the ERA5-Land reanalysis dataset. The results revealed the following:

(1)

The SPI, SSMI, and SRI exhibited a high degree of spatial and temporal consistency in their trends with an insignificant decreasing trend dominating in most parts of the upper region and an insignificant upward trend dominating in most parts of the lower region. Moreover, the SRI displayed the most significant variation, with most regions showing a significant increasing/decreasing trend.

(2)

Both linear and nonlinear methods exhibited strong temporal and spatial consistency under PMAD and PMHD, with linear relationships being stronger than nonlinear ones. More than 80% and 70% of the study area showed identical DRTs for the two methods for PMAD and PMHD, respectively.

(3)

The DRTs of PMAD and PMHD were around 1–2 months and 3–5 months, respectively. Significant differences existed in the DRT between the dry season and the rainy season. For agricultural drought, the DRT was 1 month in the dry season and 1–2 months in the rainy season. Regarding hydrological drought, the DRT was 1–3 months in the dry season and 3–5 months in the rainy season.

(4)

Divergent spatial patterns of the proportion of DRT were observed between PMAD and PMHD. The upper sub-regions with a larger proportion of areas showed a longer DRT of PMAD but a shorter DRT of PMHD, while the lower sub-regions with a larger proportion of areas showed a shorter DRT of PMAD but a longer DRT of PMHD.

(5)

Significant statistical correlations between meteorological drought and agricultural drought and between meteorological drought and hydrological drought were observed in specific periods for each sub-region. Significant coherence was exhibited between SPI-m (m represents the DRT under PMAD in each sub-region) and SSMI and between SPI-n (n represents the DRT under PMHD in each sub-region) and SRI in most parts of each sub-region. This suggests that meteorological drought is the key driver of agricultural and hydrological drought.

(6)

Hydrometeorological factors and environmental characteristics collectively influenced the DRT. The hydrometeorological factors contributed the most to DRT, followed by terrain factors, and the land cover types contributed the least. Specifically, the areas with increased precipitation, soil moisture, or runoff had a shorter DRT, and the regions with higher elevations and slopes tended to exhibit a longer DRT of PMAD and a shorter DRT of PMHD, while the areas with a low proportion of cropland and high proportion of forest tended to display shorter DRTs.