Investigating the vegetation's temporal–spatial response to meteorological and hydrogeological drought in drylands

As the most constraining environmental factor of vegetation development in arid areas, soil moisture is mainly supplied by precipitation and groundwater resources. Considering the varying sensitivity of different plant communities to precipitation and groundwater‐induced water shortages, the communities' management requires the study of the effects of both meteorological and hydrogeological droughts on vegetation. Therefore, this study sought to model the effects of meteorological and hydrogeological droughts on vegetation indices obtained from MODIS satellite images in the Sirjan plain from 2000 to 2019. To this end, vegetation communities were first identified and separated based on extensive field operations, taking into account the Normalized Difference Vegetation Index (NDVI) and Vegetation Condition Index (VCI) at the plant communities' scale. Then, meteorological drought was calculated using Standardized Precipitation Index (SPI) and the hydrogeological drought was also measured by Groundwater Resource Index (GRI) via the Kriging technique. Finally, the relationship was modelled using Multivariate Linear Regression (MLR). The results revealed that SPI at a 6‐month time scale (as an important index) played a significant role in elaborating the changes in VCI in eight communities out of 18 ones at the 6‐month scale, acting as a strong and reliable estimator. Moreover, 61.6% of VCI changes in Artemisia sieberi‐Seidlitzia rosmarinus community were explained with GRI and SPI at a 6‐month time scale, indicating the dependence on groundwater and precipitation. Nonetheless, some communities (i.e. six cases) were unresponsive. These findings draw our attention to the importance of considering the special reaction of each plant community in the management.


| INTRODUCTION
As a complex natural disaster, drought affects a wide range of phenomena (Rimkus et al., 2017), bringing about diverse negative environmental and socioeconomic consequences (Huang et al., 2017;Mishra & Singh, 2010).Moreover, compared with other natural disasters, this phenomenon includes more impact time and a higher spatial range (Bagheri & Mohamadi, 2012).Therefore, the severity of the negative consequences of drought is expected to increase within the years to come, especially in mid-latitudes and arid regions (Dai, 2011).On the contrary, drought is difficult to assess and monitor due to such factors as the inseparability of its beginning and end, and the lack of a clear universal definition (Ji & Peters, 2003).However, drought indices including the meteorological, agricultural and hydrological ones can readily convey relevant information to different audiences, with meteorological indices being able to explicate climatic anomalies in terms of severity, duration and spatial range (Wilhite, 2000).McKee et al. (1993) developed the meteorological index of SPI, obtained by correlating the precipitation statistics and the gamma probability distribution over a specific time and place.Based on the purpose of the study, SPI can be measured and assessed multidimensionally over time, enabling the monitoring of the effect of drought on different ecosystem components (Caparrini & Manzella, 2009).Furthermore, the simplicity of the SPI has led to its widespread use (Stagge et al., 2015).As a reliable index for assessing the onset of drought, SPI is an important tool in assessing the moisture conditions and realizing the responses in the ecosystem (Ji & Peters, 2003).However, using merely one index to assess the impact of drought on complex ecological systems is not recommended (Steinemann & Cavalcanti, 2006).While groundwater resources are one of the most important elements of the hydrological cycle (Nielsen & Ball, 2015), their effect on vegetation is usually less obvious and more difficult to assess (Zhu et al., 2020).On the contrary, groundwater is mainly fed by precipitation, exerting a great influence on the temporal and spatial distribution of the soil moisture on which vegetation depends (Dai et al., 2014;Naumburg et al., 2005).The relationship between NDVI and the groundwater's depth has already been studied on local scales (Doble, 2006;Jin et al., 2014;Zhu et al., 2015).However, it is very difficult to study the function of groundwater (feeding and discharging), as it requires a threshold to detect the phenomenon of drought in groundwater (Elbeih, 2015;Jha & Chowdary, 2007).Moreover, the groundwater's feeding, level and discharge decrease with the occurrence of drought, respectively (Roshun & Habibnejad Roshan, 2018).Applying hydrogeological drought indices such as GRI can help examine the changes in this variable (Mohammadi et al., 2018;Sharifi et al., 2020).GRI was introduced and assessed by Mendicino et al. (2008) as a drought index.Standardizing the piezometric wells' data concerning the groundwater level changes based on the Z number, hydrological GRI acts similarly to the SPI (Khosravi et al., 2016).
Vegetation is considered as a sensitive index in studying global climate change, especially in dry areas (Kutiel et al., 2004;Li et al., 2003) where, on the one hand, there is a great dependence on groundwater and the changes in the surface of the groundwater aquifer lead to temporal and spatial changes in vegetation (Dai et al., 2014), and, on the other hand, there are considerable changes in important seasonal and annual climatic indices such as precipitation (Sonfack et al., 2013).Therefore, meteorological and hydrogeological drought indices are useful tools for showing changes in the vegetation's water resources.However, the vegetation response to these characteristics varies in different vegetation communities (Rousta et al., 2020).The vegetation's response pattern to these indices may vary due to different root systems, location and meristem density of vegetation communities' dominant species.Therefore, it appears necessary to consider the vegetation communities as the working basis and a homogenous area.NDVI (as an index for green biomass), leaf area index, production, greenness and the ecosystem's photosynthetic capacity have already been attended to by various researchers (Zhang et al., 2013;Zhao et al., 2018).VCI (which is related to NDVI) has also been considered by scholars due to its ability to exclude the effect of local geographical factors on vegetation (Quiring & Ganesh, 2010).
While the relationship between vegetation indices, meteorological drought (Jain et al., 2010;Li et al., 2015;Xu et al., 2018;Zhang et al., 2017;Rimkus et al., 2017;Li et al., 2018;Zhao et al., 2018;Barker et al., 2016;Niu et al., 2019;Rousta et al., 2020;Zhong et al., 2020;Kaur et al., 2022) and hydrogeological indices (Doble, 2006;Jin et al., 2007Jin et al., , 2014) ) has already been studied separately, the concurrent relationship between vegetation indices, meteorological drought and groundwater has been underresearched, especially at the scale of plant communities (e.g.Jin et al., 2016), which is one of the objectives of this study.Having sufficient knowledge about the response of each plant community to drought indices is crucially necessary.It is a major area of interest within the field of the crisis management at both pre-crisis and post-crisis stages.As homogeneous vegetation, each plant community has its own response pattern to droughts due to limited available water resources (Bagheri et al., 2021), which can be examined through the specular reflection of the vegetation indices obtained from satellite images.
Up to now, far too little attention has been paid to investigate concurrent relation of meteorological and hydrogeological droughts with vegetation in the plant community scale.This indicates a need to understand the simultaneous roles of main limiting variables that lead to changes in vegetation of drylands as dependent one which is attractive for ecologist and executive managers in terms of applying these reactions in the management and decisions at the drought crisis period.Therefore, a concurrent consideration of meteorological and hydrogeological droughts as independent ones using the stepwise Multivariate Linear Regression (MLR) models was applied to predict the changes in vegetation indices as the responds in plant community scale.It should be noted that the concurrent consideration of many independent variables to explore influencing factors is one of the most important advantages of the MLR models in environmental sciences (e.g.Cho & Lee, 2018).
As the main influential factors on vegetation in arid regions, precipitation and groundwater are the primary sources of soil moisture (Pereira et al., 2002).Considering the different sensitivity of plant communities to the precipitation and groundwater-induced water scarcity caused by droughts (Páscoa et al., 2020;Rousta et al., 2020), the study of the combined effects of meteorological and hydrogeological drought on vegetation becomes considerably important for the management of plant communities.Therefore, this study was carried out to model the influence of meteorological and hydrogeological drought on vegetation indices of plant communities obtained from MODIS satellite images over a statistical period of 20 years (2000-2019).

| Study area
Study area (i.e.Sirjan plain with the plain's code of 4,419) covering 3,718.37 km 2 is located in the central plateau basin in Kerman province of Iran.It has been extended between east longitude 55°15 and 56°10′ and north latitude 28°56′ and 29°50′ (Figure 1).It is considered one of the watersheds of the Abarqo-Sirjan Kavir (as a closed basin in the central plateau of Iran), situating in the vicinity of the Kafa Namak playa.Annually, mean precipitation equals 113.5 mm and yearly average temperature is 18.17 Celsius.The vegetation is mainly included by shrubs and bushy trees.

| Research method
In order to study the response of vegetation to meteorological and hydrogeological droughts, the research process is demonstrated in the following flowchart (Figure 2).
Since the modelling of vegetation index response in each plant community to drought and hydrogeological indices was aimed in this study, first, plant communities were identified based on extensive field operations matched with June month having the maximum vegetation.The different plant communities were named according to the dominated species.In this regard, the name of the first dominant species was determined based on the highest percentage of canopy cover, and then, the next dominant species were added in naming the plant community if they had more than half of the canopy cover (%) of the dominant species of the previous stage.
Vegetation indices were then extracted MODIS satellite images from NASA's site at the time of year associated with the maximum vegetation (i.e.June).The NDVI as one of the most widely used methods for investigating vegetation (Zhao et al., 2018) and the VCI as an index eliminating the effect of local geographical factors on NDVI index (Quiring & Ganesh, 2010) were calculated based on the following formulas: where NIR is related to near-infrared band and RED is belonged to red band and also NDVI min and NDVI max are related to minimum and maximum NDVI levels during the 20-year study period in each plant community, respectively.
The NDVI and VCI values of plant communities were obtained via limiting the maps based on communities' boundaries using Arc-GIS 10.4.1 software environment.In this regard, a practical advantage of using the VCI is that it excludes the effect of local geographical factors on vegetation.
After determining the study period to assess the meteorological drought in the study area, the representative station (i.e.Sirjan synoptic station) was selected.In this regard, the Standardized Precipitation Index (SPI) was applied in periods of 1, 3, 6, 9, 12, 16 and 24 months, as well as annually using the MATLAB software environment on monthly data of precipitation from 2000 to 2019.It should be noted that gamma probability density function is fitted to a given frequency distribution of precipitation series in a specific station (Liu et al., 2021) based on McKee's theory in 1993.Accordingly, SPI can then be calculated based on the following formulas: where α is a parameter about shape, β is a parameter about scale and x equals the amount of precipitation.It should be noted that the best values of α and β are estimated using the maximum likelihood method and the gamma function is displayed as: wherein n means the number of precipitation series.
The cumulative probability for a given month can thus be acquired using the following equation: the SPI can be calculated as follows: wherein x is the amount of precipitation and G(x) is Γ function-related precipitation probability distribution, S is the positive and negative coefficient of cumulative probability distribution, when G(x) > 0.5, S = 1 and when G(x) ≤ 0.5, S = −1.c0 = 2.5155, c1 = 0.8028, c2 = 0.0103, d1 = 1.4327, d2 = 0.1892, d3 = 0.0013.Finally, the SPI categories are presented in Table 1.
For the purpose of hydrogeological drought analysis, the first step in this process was to prepare the point data of groundwater level of piezometric wells in the plain surface from Kerman Regional Water Company.The second step set out to the temporal and spatial zoning of groundwater level data in Gs + software package using Kriging geostatistical technique as a strong estimator (Bagheri & Mohamadi, 2012) with the best variogram model fitted data based on the following formula: wherein: γ(h): The value of the variogram for a pair of points spaced h apart.n(h): The number of pairs of samples used for a given distance such as h.Z(x): The variable observed at point x.Z(x + h): The observed value of a variable that is at a distance of h from x.
The third step sought to determine the groundwater level of each plant community via limiting the spatial zonation maps of the groundwater level of the plain to the map of plant community boundaries in the Arc-GIS software package.Furthermore, the nonparametric Mann-Kendall test (Kendall, 1975) associated with Theil-Sen's estimator having confidence level lines (Cao et al., 2014) was carried out on the plain's groundwater table time series for assessment of trend.
Finally, the groundwater level index of each plant community was normalized by applying the Groundwater Resource Index (GRI) according to the below formula:  et al., 1993;Mendicino et al., 2008).
where X ij is the groundwater level for the i unit (in this study the plant community) and j observation, X m belongs to the mean groundwater level and SD relates to the standard deviation.After limiting zonation maps of groundwater table from 2000 to 2019 with plant communities including 18 zones (Figure 3) via zonal statistics in Arc-GIS 10.4.1 software package, all statistical characteristics were extracted to calculate the GRI of each plant community under a normal statistical distribution.The criteria for evaluating the GRI of plant communities were based on the changes in the GRI classes in a certain time.In this regard, standard classes was presented in Table 1.
To end the process, modelling of relationship between drought indices (i.e.SPI and GRI) as stimuli and vegetation index (i.e.VCI) as a dependent variable at the plant community scale was performed adopting multivariate regression method with stepwise procedure in SPSS software (version 24).The validation of the modelling was carried out via R 2 criteria, which can be calculated from the following equation: where in Q obs, i is the observable value of the VCI, Q est, i relates to the estimated value of the VCI, Q obs,i belongs to average observational value of VCI, Q est, i is the average estimated value of the VCI and i is the time step.

| Vegetation indices
The results of the field operation for determining plant communities' map (i.e.typology) are shown in Figure 4, indicating the recognition of 18 plant communities in the As shown in Figure 5, most vegetation changes in vegetation communities were not affected by time.Accordingly, changes in NDVI of Alhagi camelorum and Cornulaca monacantha-Salsola orientalis communities were the most influential from time (Figure 5a, f) and conversely change in NDVI of Ephedra pachyclada-Noaea mucronate one was the least effective (Figure 5o).
The results of the VCI analysis are displayed in Figure 6.It can be seen from VCI in Figure 6 that Seidlitzia rosmarinus-Artemisia sieberi community changed dramatically between years and conversely Seidlitzia rozmarinus-Salsola tomentosa one had the least between-years variations of VCI with amount of 0.554.This difference might be attributed to the stability of second-dominated plant of these communities over droughts and other environmental influencing factors.

| Meteorological drought
The results of meteorological drought index evaluation in 1-, 3-, 6-, 9-, 12-, 16-and 24-month time bases are presented in Figure 7 (Figure 7a, g).As can be seen from the Figure (below), the most severe dry month of 1-month scale occurred in the first month of 2013 (Figure 7a) and also the wettest 3-month one happened in 2013 (Figure 7b).On a 6-month scale, the wettest case occurred in 2013 (Figure 7c).At the 9-and 12-month scales, the driest case is observed in 2008 (Figure 7d, e).Finally, on a scale of 16 and 24 months, the driest cases were found in 2001 (Figure 7f, g).
In addition to the above, the SPI index was plotted as an annual SPI based on annual rainfall changes (Figure 7h).
What stands out in figure is the existence of highseverity events in droughts.In this regard, although the SPI of 2004SPI of , 2011SPI of , 2013SPI of , 2015SPI of , 2017SPI of and 2019 was in the middle-wet class, but this index for the year of 2001 became in the severe-drought class and for the years of 2000, 2001, 2008, 2010 and 2016 gained in the moderatedrought class.

| Hydrogeological drought
The results of groundwater level index assessment in the plain level by zoning method with geostatistical estimator in the studied years (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018)(2019) are displayed as maps in Figure 8.According to these results, the groundwater level from 2000 to 2019 in different parts of the plain had irregular changes.Closer inspection of the maps in Figure 8 shows that the groundwater level in the northern parts of the plain is higher than the rest of the plain in all years.Also, the western and southern parts of the plain have the lowest groundwater level in all years.The results of trend assessment for plain's groundwater table time series using the Mann-Kendall test (Figure 9) showed a significant negative trend (i.e.p < 0.01) for the studied plain from 2000 to 2019 due to the fact that the amount of z exceeded of ±1.96.Accordingly, Theil-Sen's test revealed that 1 m of water table drop (on average) can be observed in the studied plain every year, indicating abrupt changes in some years including 2000, 2001, 2002, 2006, 2009 and 2019.Irregular temporal and spatial changes in water table evaluation created a necessity to calculate the GRI for different plant communities for each year.Therefore, the results of GRI as an index of hydrological drought for different plant communities are set out in Figure 10.The most interesting aspect of this figure is dramatical changes in GRI for different plant communities.In terms of allocating the most frequently negative GRI, Ephedra pachyclada-Noaea mucronate having 12 years of 20 ones was recognized as the first (Figure 10o), followed by Artemisia sieberi-Seidlitzia rosmarinus and Cornulaca monacantha-Salsola orientalis having negative GRI as 11 years of 20 ones are in the second grade (Figure 10b, f).

| Regression modelling
The results of modelling the vegetation response to meteorological and hydrological drought using MLR are presented in Table 2 with the validation results in Figure 11.As can be seen from the table (below), the meteorological drought index played an important role in explaining the VCI changes in all studied plant communities, while the  effect of the GRI hydrogeological drought index was found in explaining the VCI changes in only one community (i.e.Artemisia sieberi-Seidlitzia rosmarinus).In this regard, a strong and positive relationship was established between only SPI drought index in 6-month time scale and VCI index in some plant communities.In addition, the 9-month SPI drought index along with the 6-month SPI one explained 83.8% of the VCI changes in Cornulaca monacantha-Launaea procumbens community.The SPI in 9-month and annual time scales associated with SPI in 6-month time scale explained 88.1 per cent of VCI changes in Astragalus spachianus plant community and SPI in 6month time scale in addition to the annual SPI explained 71.4 per cent of the VCI changes in plant community of Salsola arbusculiformis-Zygophyllum eurypterum.

| DISCUSSION
The results of vegetation response to meteorological and hydrogeological drought revealed that only 6-month SPI as the most important indicator played a major role in explaining the VCI changes in eight plant communities including Artemisia sieberi-Zygophyllum eurypterum, Cornulaca monacantha-Launaea acanthodes, Cornulaca monacantha-Salsola orientalis, Dendrostellera lessertii-Noaea minuta, Halocnemum strobilaceum-Seidlitzia rosmarinus, Dendrostellera lessertii-Noaea mucronata, Ephedra intermedia-Cousinia multilobe and Seidlitzia rosmarinus-Salsola dendroides.So that this index alone was a strong and reliable estimator for estimating the VCI of the mentioned plant communities.Additionally, 83.8% of the VCI changes in the plant community of Cornulaca monacantha-Launaea procumbens were attributed to the SPI in 6-month time scale associated with the SPI in 9-month time scale.However, SPI in 6-month and annual time scale along with SPI in 9-month time scale as explanatory factors explained 88.1 per cent of VCI changes in Astragalus spachianus plant community and SPI in 6-month time scale along with annual SPI drought index caused 71.4 per cent of the VCI changes in Salsola arbusculiformis-Zygophyllum eurypterum.
These results indicate the response of vegetation of some plant communities to the meteorological drought index.This result was in accord with the results of Quiring and Ganesh (2010) indicating a strong correlation between the 6-and 9-month SPI with the VCI.In a study of the relationship between SPI at time scales of 1, 3, 6, 9 and 12 months with NDVI in three regions of India, it was found that in each region, depending on the growth rate of vegetation, the optimal relationship occurred at different time scales (Jain et al., 2010).The best explanatory factor of SPI at different time scales for NDVI was reported SPI in 6-month scale (Caccamo et al., 2011).Farrokhzadeh et al. (2017) claimed that the correlation coefficients between SPI with NDVI and EVI are 0.78 and 0.75, respectively.In the study of vegetation response to meteorological drought in different biomes in northern China (i.e.vegetation of temperate steppes, temperate desert steppes, shrublands and arid forests), Xu et al. (2018) showed a significant relationship between vegetation indices and meteorological drought.In the same vein, in a study investigating vegetation reaction to SPI, Fang et al. (2019) reported that vegetation deterioration probabilities under meteorological drought reached their peak (39.7%) in summer, indicating the highest vegetation vulnerability to drought stress.In a study conducted by Zhong et al. (2020), it was shown that VCI with SPI had a significant lag (i.e. 6 months) correlation in the upstream Heihe River Basin, Northwest China.The relation between vegetation indices (i.e.NDVI and VCI) and meteorological drought index (i.e.SPI) in the majority of plant communities in the present study are in accord with recent studies including Quiring and Ganesh (2010), Jain et al. (2010), Caccamo et al. (2011), Farrokhzadeh et al. (2017), Xu et al. (2018), Fang et al., (2019), Zhong et al. (2020) and Wei et al. (2022).Owing to the possibility of extracting drought time series through satellite images and establishing a strong  The results of this study showed significant differences in the GRI in different plant communities    during different years.Because in the study area, where spring rains played a significant role in the growth of maximum plant communities, the alignment between changes in the hydrogeological drought index (i.e.GRI) and VCI was observed only in phreatophyte plant communities.In this regard, 61.6 per cent of VCI changes in Artemisia sieberi-Seidlitzia rosmarinus community were explained by GRI hydrological drought index along with SPI meteorological drought index in 6-month time scale, indicating the dependence of this plant community on groundwater in addition to precipitation.This valuable result of the present study designed to determine the response of communities to different types of droughts (i.e.meteorological and hydrogeological ones) may help us to enhance our understanding of explanatory factors of the changes in vegetation via VCI in June when vegetation growth is maximum.However, in order to investigate the relationship between only hydrogeological drought and vegetation, it needs to assessment in dry seasons (i.e.mid-summer), when there is no precipitation and soil moisture plays a limiting role in the growth of communities.Furthermore, in-depth investigations between the hydrogeological drought index (i.e.GRI) in plant communities and other characteristics (i.e.reproductive indices) are recommended to examine degree of plant community's reliance on groundwater reserves.Since our goal was to explore the strength of meteorological and hydrogeological drought-based model for each plant community (as study unit) in the region, validating the reliable models (i.e. with an R 2 more than 0.6) is recommended in the neighbouring areas having similar communities.On the contrary, some plant communities (i.e. six cases) including Alhagi camelorum, Seidlitzia rosmarinus-Artemisia sieberi, Artemisia sieberi-Zygophyllum atriploides, Ephedra pachyclada-Noaea mucronate, Haloxylon ammodendron and Seidlitzia rozmarinus-Salsola tomentosa had insignificant and in one case, unreliable relationships with meteorological and hydrological drought in the MLR models (i.e.having R 2 less than 60%), indicating the stability and immutability of these communities.Therefore, preservation of these constant plant communities is one of the recommendations of this research to the executive bodies and officials of natural resources.

| CONCLUSION
The purpose of the current study was to explore the influence of meteorological drought (SPI) and hydrogeological one (GRI) on vegetation index (VCI) in the plant community scale.Our finding revealed that (1) some plant communities (i.e. 8 cases) are changed significantly by only meteorological drought (i.e.SPI-6 month), indicating feedbacks to meteorological drought as responses.(2) One plant community (Artemisia sieberi-Seidlitzia rosmarinus) is affected significantly by meteorological and hydrogeological drought, indicating sensitive case to both types of droughts.(3) Some plant communities (i.e. six cases) had insignificant relationships with meteorological and hydrological drought, indicating apathetic and unresponsive communities.
These results draw our attention to the importance of applying the special reaction of each plant community as a manageable tool by human over threatening factors including meteorological and hydrogeological droughts in the natural ecosystems of drylands.Since the meteorological and hydrogeological droughts induced by climate change can leave deteriorating effects on the vegetation and greenness, implementing the conservation, re-habitation and management projects in natural ecosystems by the executive managers and decision-makers are then recommended in the present and future years on the basis of exploring and applying the reactions to have more native greenness in drought crisis period.Further studies require to better understand the special behaviour of plant communities in drylands against the meteorological and hydrogeological droughts across the desert biome as the most vulnerable zone to climate change and its consequences (i.e.drought).

F
Location of study area in Kerman province of Iran (a) associated with its piezometers and meteorological site as (Sirjan synoptic station) (b).
Zones used for implementing zonal statistics in the current study.studied area.The findings of the NDVI assessment for all recognized plant communities with different vegetative forms including forb, plant and shrub in the study years (2000-2019) are shown in Figure 5. Looking at Figure 5, it is apparent that the highest amount of NDVI was observed in the plant community of Dendrostellera lessertii-Noaea minutes, which indicates the high greenness of this plant community and the lowest amount of NDVI was related to Haloxylon ammodendron with amount of 0.070.
Groundwater level maps provided by Kriging estimator in the Sirjan plain from 2000 to 2019.

F
Trend assessment for plain's groundwater table during 20 years (2000-2019) using the Mann-Kendall test and Theil-Sen's estimator.
Flowchart illustrating the processes of conducting the research.
T A B L E 1 GRI (or SPI) classification for drought assessment (McKee