Deep Learning for Monitoring Agricultural Drought in South Asia Using Remote Sensing Data

: Drought, a climate-related disaster impacting a variety of sectors, poses challenges for millions of people in South Asia. Accurate and complete drought information with a proper monitoring system is very important in revealing the complex nature of drought and its associated factors. In this regard, deep learning is a very promising approach for delineating the non-linear characteristics of drought factors. Therefore, this study aims to monitor drought by employing a deep learning approach with remote sensing data over South Asia from 2001–2016. We considered the precipitation, vegetation, and soil factors for the deep forwarded neural network (DFNN) as model input parameters. The study evaluated agricultural drought using the soil moisture deﬁcit index (SMDI) as a response variable during three crop phenology stages. For a better comparison of deep learning model performance, we adopted two machine learning models, distributed random forest (DRF) and gradient boosting machine (GBM). Results show that the DFNN model outperformed the other two models for SMDI prediction. Furthermore, the results indicated that DFNN captured the drought pattern with high spatial variability across three penology stages. Additionally, the DFNN model showed good stability with its cross-validated data in the training phase, and the estimated SMDI had high correlation coefﬁcient R 2 ranges from 0.57~0.90, 0.52~0.94, and 0.49~0.82 during the start of the season (SOS), length of the season (LOS), and end of the season (EOS) respectively. The comparison between inter-annual variability of estimated SMDI and in-situ SPEI (standardized precipitation evapotranspiration index) showed that the estimated SMDI was almost similar to in-situ SPEI. The DFNN model provides comprehensive drought information by producing a consistent spatial distribution of SMDI which establishes the applicability of the DFNN model for drought monitoring.


Introduction
Drought is regarded as a common and consistently occurring weather related phenomenon that has a severe impact on human society and ecosystem [1][2][3]. This leastunderstood natural phenomenon is very challenging to detect, with the frequency and intensity of events varying considerably due to frequent global climate change [4]. Though a drought generally starts with a precipitation shortage, depending on its mechanism as support vector regression (SVR), random forest (RF), boosted regression tree (BRT), bias-corrected random forest (BRF), gradient boosting machine (GBM), and classification and regression tree (CART). For example, Zhang et al. [30] used the XGBoost model for meteorological drought prediction. The results suggested that XGBoost had high predictive skill compared to the lag distributed non-linear model (DLNM). In another study, Chiang and Tsai [31] found SVM (Support Vector Machine) superior in predicting hydrological drought compared to the traditional model. In some cases, however, the performance of deep learning in drought studies exceeds other machine learning models [28]. However, very limited studies found in the literature used deep neural nets for drought monitoring and forecasting. Lee et al. [32] used a deep learning model for soil moisture estimation from satellite remote sensing data to monitor drought. The model showed high stability during cross-validation and high correlation with in-situ observation. Zhang et al. [33] reported that soil moisture estimated by the deep learning model from remote sensing data was able to capture the complex relationships of in-situ soil moisture. Although Agana and Homaifar [34] used a deep belief network (DBN) for long-term drought prediction, they only considered the lagged value of a standardized streamflow index (SSI) as input parameters for SSI-12 and SSI-24 prediction. Though the drought monitoring model constructed by Shen et al. [25] using a deep learning approach, however, in this study we proposed an agricultural drought predicting system using a deep learning model taking into account precipitation, soil, and vegetation as explanatory variables to monitor drought in South Asia. Moreover, the variability of agricultural drought patterns was investigated in relation to the crop phenology stage, which was not reported before in the South Asia region. Two machine learning models such as distributed random forest (DRF) and gradient boosting machine (GBM) were used to compare performance with the deep learning model. As a result, the main goal of this study is to develop an integrated drought-monitoring system considering precipitation, soil, and vegetation factors using a deep learning approach. Therefore, the present study aims to: (i) use a deep learning model for monitoring agricultural drought from 2001-2016; (ii) characterize spatiotemporal variation of agricultural drought during three crop phenology stages; and (iii) examine the performance of deep learning for predicting the SMDI.

Study Area
The study area is the southern region of Asia, covering about 5.2 million Km 2 geographically positioned at 5 • -40 • N and 60 • -100 • E [35]. The region encompasses eight countries (India, Nepal, Bhutan, Bangladesh, Pakistan, Afghanistan, Sri Lanka, and the Maldives) with about 1.836 billion population [36]. The area is defined by the Indian Ocean and topographically dominated by the Indian Plate. South Asia is covered by a heterogeneous land surface area ( Figure 1). The diverse climatic conditions that prevail in South Asia vary from region to region. For example, a hot summer exists in the southern parts, though heavy rainfall occurs during the monsoon period. The indo-geographic plain of the northern side is also hot in summer and cool in winter, while snowfall occurs in Himalayans regions in the mountainous north. In general, continental climate is observed in north India and Pakistan. India in the south and Sri Lanka in the southwest have an equatorial climate. Bangladesh's climate is largely characterized by a cool winter and hot, humid summer, while Afghanistan is dominated by dry areas [14]. Seasonal rainfall is primarily driven by monsoon wind, which accounts for 80% of the annual precipitation in most parts of South Asia. Recently, climate change has increased land surface temperatures, and erratic rainfall makes South Asia more prone to drought [37].

Data
In this study, data sets from remote sensing and ground observation (weather stations) were used to calculate the drought indices for deep learning model input parameters. The specification of these data sets is presented in Table 1. is a multispectral medium/high resolution sensor consisting of terra and aqua satellites with a wide range of applications in the field of earth and environmental science [38]. MODIS provides valuable information by detecting electromagnetic energy in a wide spectral range to study the earth's ecological, meteorological and hydrological condition [39]. In this present study, we used MODIS vegetation indices product (MOD13A1) with a spatial resolution of 500m and land surface temperature (LST) product (MOD11A2) with a 1km spatial resolution for VCI and TCI (temperature condition index) drought indices calculation. MOD13A1 products can be used for global vegetation change detection and monitoring vegetation

Data
In this study, data sets from remote sensing and ground observation (weather stations) were used to calculate the drought indices for deep learning model input parameters. The specification of these data sets is presented in Table 1.  Imaging Spectroradiometer) is a multispectral medium/ high resolution sensor consisting of terra and aqua satellites with a wide range of applications in the field of earth and environmental science [38]. MODIS provides valuable information by detecting electromagnetic energy in a wide spectral range to study the earth's ecological, meteorological and hydrological condition [39]. In this present study, we used MODIS vegetation indices product (MOD13A1) with a spatial resolution of 500m and land surface temperature (LST) product (MOD11A2) with a 1km spatial resolution for VCI and TCI (temperature condition index) drought indices calculation. MOD13A1 products can be used for global vegetation change detection and monitoring vegetation biophysical interaction with photosynthetic activity [40]. MOD11A2 product consists of day and night time LST calculated within an 8 day period by combining surface atmospheric Remote Sens. 2021, 13, 1715 5 of 28 interaction and energy fluxes considering both ground and atmosphere [41]. A number of 10 tiles (h22v05-h26v05, h23v06-h26v06m, h24v07-h26v07, and h25v08-h26v08) were downloaded from the LAADS (Atmosphere Archive and Distribution System) website ( https://ladsweb.modaps.eosdis.nasa.gov/) (accessed on 27 September 2020) overing the study area from 2001-2016.
2.2.2. GLDAS-NOAH Soil Moisture, Evapotranspiration (ET), Potential Evapotranspiration (PET) Data GLDAS (Global Land Data Assimilation System), a global high-resolution terrestrial modelling system, consists of Noah, CLM (community land model), VIC (variable infiltration capacity) and land surface model (LSM), and provides land and water energy fluxes and also different soil properties based on satellite and ground observation for ecosystem modelling [42][43][44]. We used the GLDAS NOAH (GLDAS_NOAH025_M_2.1) soil moisture, ET, and PET product for this study with a spatial resolution of 00.25 • × 0.25 • , downloaded from the Land Data Assimilation System (LDAS) website (https: //ldas.gsfc.nasa.gov/gldas) (accessed on 27 September 2020). Monthly soil moisture products of the topsoil layer at 10 cm depth were used for the SMDI (Soil Moisture Deficit Index) drought index calculations for this study. The EDI (Evaporative Drought Index) was calculated using ET and PET products to evaluate surface dryness.

CHIRPS Data
CHIRPS (Climate Hazards Group InfraRed Precipitation with Station data) is a relatively high-resolution quasi-global long term (1981-present) satellite precipitation product developed by the Climate Hazards Center for monitoring drought and conducting precipitation trend analysis [45]. CHIRPS is a satellite-estimates product blended with gauge observation from GHCN (Global Historical Climate Network) and GSOD (Global Summary of the Data set) data sources [46]. The CHIRPS data was downloaded from the Climate Hazards Center's website (https://www.chc.ucsb.edu/data/chirps) (accessed on 27 September 2020) for the 2001-2016 period. Monthly precipitation data, which is available at 0.05 • spatial resolution, was used for PCI (Precipitation Condition Index), PAI (Precipitation Anomaly Index), and SPI (Standardized Precipitation Index) calculation as satisfactory performance was found with this data set in many studies [45,47], due to consistency and sufficient length which make the data more reliable for climate variability analysis [48,49].

GIMMS-NDVI Data
In this study, the NDVI data set was derived from reflectance observed by the AVHRR (Advanced Very High Resolution Radiometer) as a part of the Global Inventory Modeling and Mapping Studies (GIMMS) project obtained from https://ecocast.arc.nasa.gov/data/ pub/gimms/3g.v1/ (accessed on 27 September 2020) with a data period of 2001 to 2015. This data set has 1/12 • × 1/12 • spatial resolution and consists of the maximum NDVI value for each 15-day interval. This data set was used to generate phenology metrics for South Asia.

Meteorological Station Data
We used ground observation data, such as precipitation and temperature, from all 780 available weather stations over the study area during 2001-2016, collected from the department of meteorology for each country. In this study, we calculated SPEI following the Hargreaves [50] method from monthly precipitation and temperature data using the SPEI package in R software. First, daily precipitation and temperature were collected, then aggregated to monthly intervals for SPEI calculation. The scales of three, six, and twelve months (SPEI-3, SPEI-6, and SPEI-12) were taken into account for the SPEI calculation used to compare the results estimates from the DFNN model.

Model Input Parameter
Drought is a relatively complex event, unlike other natural disasters, due to its slow development over a prolonged period. Drought characteristics rely on a given area's climate characteristics and are substantially different across a given region; for example, a region with an extended dry period may lead to drought. Considering its evolution process, however, drought largely depends on precipitation, ultimately affecting soil moisture. As a result, drought affects the vegetation on the ground. Deficit precipitation for a long time creates a shortage of surface and sub-surface water supply [51][52][53]. Given these complex mechanisms, drought monitoring methods need to improve, which depends largely on monitoring precipitation, soil and vegetation factors [54,55]. As a result, we adopted precipitation, soil, and vegetation factors as input parameters for our deep learning model to monitor drought in the study area. Thus, the twelve variables used for model input parameters are presented in Table 2 with detailed descriptions. × 100 (where P is the precipitation in a month, and P max and P min are maximum and minimum precipitation in a month) [30] PAI P ai = P i −p p × 100 (where P i is the precipitation in a month, and P is the average annual precipitation in a month) [25] SPI-3 SPI-6 SPI-12 (k represents precipitation probability function; c 0 c 1 , c 2 , c 3 , d 1 , d 2 , and d 3 are constant) [56,57] SPEI-3 SPEI-6 SPEI-12 w − c 0 +c 1 w+c 2 w 2 1+d 1 w+d 2 w 2 +d 3 w 3 (w is defined as climatic water balance calculated based on the difference between precipitation and reference evapotranspiration; and c 0 c 1 , c 2 , c 3 , d 1 , d 2 , and d 3 are constant) [ where SMDI j−1 represents the SMDI for first month and SM j denotes soil moisture anomaly of month j) Agricultural drought arises as a consequence of meteorological drought due to insufficient precipitation, which in turn results in low water content in soil [33]. Agricultural drought was measured by SMDI drought index in this study, which was considered as observed SMDI. Moreover, ET is a very important factor that represents the status of soil moisture availability. Thus, surface dryness in response to the soil moisture was evaluated through the EDI drought index based on ET and PET representing the soil factors. The precipitation factor greatly influences meteorological drought which was calculated using PAI and PCI drought index in this study [30]. Besides these indexes, SPEI [58] and SPI [64] were calculated on a 3-month, 6-month, and 12-month time scale, which are able to reflect the dry and wet condition of a particular area. SPI and SPEI were confirmed as preferable meteorological drought indices by several researchers in their studies [65][66][67] due to their simplicity and flexibility for calculation at different time scales. Furthermore, crop growth is hindered by high surface temperature, which was measured in this study using the temperature condition index (TCI) [68]. Vegetation shows its response to a state of low precipitation and deficit soil moisture [69]. This situation was measured using the VCI and VHI drought indexes. Therefore, TCI, VCI, and VHI were used to represent the vegetation factor for drought characterization. We used SMDI as a response variable for the DFNN model, considering other factors (precipitation, soil and vegetation) as predictor variables. The detailed description and method of SMDI can be found in Narasimhan and Srinivasan [26]. SMDI calculation values lie between −4 to 4, which correspond to the estimates of the value from SPI. Hence, we classify drought events according to Carrão et al. [70], who followed Mckee et al. [56] drought classification ( Table 3). The detailed methodological flow chart of the study is presented in Figure 2. First, the drought indices (which represent both response and predicted variables) were calculated based on extracted phenology metrics. Then the drought monitoring model was constructed using DFNN architecture and the output of the model (predicted SMDI) was compared with two machine learning models. We validated the model with cross validation data in terms of R 2 , RMSE (root mean square error), MAE (mean absolute error) and residual deviance. The output of the model was evaluated with in-situ SPEI.

Phenology Extraction
Phenology metrics are very important indicators for ecosystem function, as delayed or advanced phenology metrics result in increased or decreased net primary production (NPP) which influences crop yield. Phenology mainly refers to the calculation of Start of Season (SOS), Length of Season (LOS), and End of Season (EOS) [71]. These phenology stages, however, are affected by varying degrees of drought, and we investigated the se-

Drought Class SMDI Range
Extreme Dry

Phenology Extraction
Phenology metrics are very important indicators for ecosystem function, as delayed or advanced phenology metrics result in increased or decreased net primary production (NPP) which influences crop yield. Phenology mainly refers to the calculation of Start of Season (SOS), Length of Season (LOS), and End of Season (EOS) [71]. These phenology stages, however, are affected by varying degrees of drought, and we investigated the severity of drought for each phenology stage. As each crop growth stage has different water requirements-for example, vegetative, flowering, and grain filling stages are more sensitive to water-we extracted the period of SOS, LOS, and EOS to quantify the water deficit during these periods. Therefore, we decided to use the same periods of each phenology stage to calculate all drought indices. The phenology metrics were calculated in R software using the greenbrown package [72]. For example, the SOS, LOS, and EOS extracted from the smoothed time series NDVI in R software varies from March to May, May to August, and September to October, respectively, from 2001-2016 over South Asia. Thus, monthly drought indices values (VCI, TCI, VHI, SMDI, EDI, PAI, PCI, SPI-3, SPI-6, SPEI-3, SPEI-6) between SOS, LOS, and EOS for each pixel were used to obtain individual years' drought status. To calculate the phenology metrics, the seasonality parameters were checked and extracted using a fast Fourier transformation of the time series NDVI data from GIMMS-NDVI (Global Inventory Modeling and Mapping Studies-NDVI). The tsgf (temporal smoothing and gap-filling) method was applied for temporal smoothing, which uses the derivative method taking into account the minimum mean annual value of each grid cell to determine SOS, LOS and EOS. Average phenology metrics extracted from the smoothed NDVI are presented in Figure 3.

Construction of Deep Learning Model Using DFNN
Deep learning is recognized as part of machine learning, capable of solving complex tasks in the form of object detection, feature selection, image classification, and decision making [73]. Deep learning varies from conventional neural networks in that it utilizes many layers and several hidden layers, whereas a traditional neural network utilizes just two to three layers [74,75]. In the present study, we used the H2O package, which was designed for DFNN as a deep learning model to use in the R software [74]. DFNN, a supervised learning model with error back-propagation, can be found in detail at

Construction of Deep Learning Model Using DFNN
Deep learning is recognized as part of machine learning, capable of solving complex tasks in the form of object detection, feature selection, image classification, and decision making [73]. Deep learning varies from conventional neural networks in that it utilizes many layers and several hidden layers, whereas a traditional neural network utilizes just two to three layers [74,75]. In the present study, we used the H2O package, which was designed for DFNN as a deep learning model to use in the R software [74]. DFNN, a supervised learning model with error back-propagation, can be found in detail at https: //docs.h2o.ai/h2o/latest-stable/h2o-docs/data-science/deep-learning.html (accessed on 27 September 2020). DFNN, also known as Feed Forwarded Neural Network (FFNN), approximates the function by learning the value of the parameters. In this way, information passes from the input nodes through hidden layers to the output layers, indicating a one-way direction with no loops in the network. Our model was a multi-layer perceptron interconnected in a feed-forward way with one input layer containing the 12 explanatory variables (predictor variables), 4 hidden layers each of 300 nodes, and an output layer (response variable) (Appendix A). A drought-monitoring model was constructed based on DFNN considering SMDI as the function of precipitation, vegetation and soil factors as presented in Equation (1). Therefore, we employed 12 drought indices as input parameters to reflect the complex interaction between precipitation, vegetation, and soil factors.
where Y = SMDI; precipitation = PAI, PCI, SPI-3, SPI-6, SPEI-3, SPEI-6; vegetation = VCI, VHI, TCI; and soil = EDI. In the present study, for each input variable a raster image of 7818 (row) × 8859 (column) dimension with a total pixel cell of 69259662 was used. Each pixel cell of a raster image represents the value of a corresponding input variable. Thus, the total size of the data is 69,259,662 for each input variable was used in the deep leaning model. First, the data set was randomly split into three parts, i.e., (i) training data (75% of the total sample), (ii) validation data (12.5% of the total sample), and (iii) test data (12.5% of the total sample). Generally, a training data set is used to fit the model with weights and biases in the case of a neural network. The model also learns from the training data set. Validation data is a sample of data used for unbiased evaluation of model fit while tuning model parameters on the training data set. Alternatively, the evaluation of the final model fit is usually done with a test data set. We applied the same data partition for 16 years (2001-2016). We calibrated the model by tuning different parameters of the model on the training data set, and the model results were validated using cross validation (a resampling procedure of the data set used to validate the model in order to estimate how the model is expected to perform) data. Further, we statistically compared the model's performance on the test data set.

Model Calibration on Training Data Set by Tuning Model Parameter
To obtain the optimal structure of a DFNN model, the number of hidden layers and nodes setup is a crucial factor, as the model's performance and accuracy largely depends on the size of the data set and model parameters. First, we used the default parameter of the H2O package to build neural networks which include the 'Rectifier' activation function with 200 neurons and 10 epoch (number of passes over the training data set to be carried out). In addition, the L1 and L2 regulation (representing the sum of square of all weights and biases in the network) used to prevent overfitting were assigned to a default value of 0. The hidden drop-out ratio that permits a large number of the models to be averaged as an ensemble was fixed to a default value of 0. Then, we examined model performance with the above-mentioned input parameters by enabling early stopping, which stops training when the model reaches a certain validation error. In order to obtain low RMSE, MAE, and mean deviance with a high R 2 value, we increased the number of neurons from 200 to 300, L1 and L2 0 to 1 × 10 −6 , hidden drop out ratio of 0 to 0.2 after each iteration, with epoch remaining at 10. The 'max_W2' (a maximum sum of squared incoming weights into anyone neuron) parameter, useful for unbound activation function, was also added with a value of 10. The 'max_W2' function introduces bias into parameter estimates, but frequently produces substantial gains in modeling as estimate variance is reduced, which enables high predictive accuracy. Moreover, the activation function was changed from 'Rectified leaner unit' to 'RectifierWithDropout' ('Rectified leaner unit' and 'RectifierWithDropout' represent non-linear activation functions) to activate the neurons more actively. We obtained minimum error for each of the performance metrics (RMSE, MAE, and deviance) when the data set irritated 10 times with 300 neurons, presented in Figure 4. We highlighted the scoring statistics considering drought severity for the year of 2012 as the same input parameters were applied for the rest of the years, hence why discussion of all years is not necessary. Further, we demonstrated the scoring history of the DFNN model at three phenology stages. Though there were some differences, the scoring trend both in training and validation data were similar. We observed, during the SOS stage, after starting the model, that the error increased sharply then decreased sharply, with some upward and downward curves as the model progressed (Figure 4a). Error became flat at epoch 10 for both data sets, demonstrating a low error rate. However, quite the opposite scenario was found at the LOS stage, which showed a sharply decreasing curve after starting the iteration. It became somewhat flat at epoch 10, though an upward and downward curve is observed between 3-8 epochs (Figure 4b). In the case of the EOS stage, the model initially showed a similar pattern to that seen in SOS, and then after the 10 th learning time, the error values were likely at a minimum compared to other learning times ( Figure 4c). After calibration of the DFNN model, we obtained R 2 for training and validation data set, as shown in Table 4. Our model produced low R 2 for training and validation data set with fewer epochs, and a gradually higher R 2 value occurred with increasing epoch numbers. The results of the model had a good agreement with R 2 values of 0.811, 0.891, and 0.793 during SOS, LOS, and EOS stages, respectively, with epoch number 10. These findings indicate that a model with a large number of epochs facilitates better learning of the training data set and, in turn, better prediction results.

Model Validation Using Cross-Validation Data
We validated the DFNN model internally using a cross-validation method that permitted us to examine how well a model learned, i.e., the stability of the network structure. A number of 5 to 10 is considered good for K-folds cross-validation, and we selected the value for K as 10. The data were randomly split for 10 fold cross-validation. All 10 cross-validation models were formed based on 90% of the training data and 10% of the validation data. The model computed validation metrics for every 10 cross-validations by scoring against the true levels of 20% validation data. As a result, to make one prediction, 10 validation predictions are merged, and in this way overall cross-validation metrics are computed. The outcome of the DFNN model cross-validation metrics is presented in Figure 5.

Machine Learning Model
Apart from deep learning, two machine learning models, DRF and GBM, were used. These two machine learning algorithms are flexible and robust for classification and regression tasks [16]. DRF is a tree-based algorithm that produces a forest of trees rather than a single tree that uses a randomized subset of candidate features by applying bootstrap sampling to make a final prediction [51]. GBM is also tree-based, using forward-learning ensemble methods that build a regression tree on all features of the data set for building a predictive model. It works by optimizing the loss function using a weak learner to make a prediction, and by applying an additive model the loss function is minimized. We used the H2O package for DRF and GBM analysis in the R software. To obtain optimum performance, we kept the same value for some important parameters such as 'number of trees', 'max. depth', 'learn_rate', 'sample_rate', fold_assignment = "Random" etc. for both the DRF and GBM models. The important parameters with their choice of values are present in Appendix B. ure 5. Figure 5 demonstrated high R 2 value ranges from 0.77 to 0.80, 0.86 to 0.90, and 0.75 to 0.80 for the periods of SOS, LOS, and EOS phonology stages, respectively. A high R 2 value indicates a low error rate which is presented as RMSE (mean value for SOS = 0.32, LOS = 0.31, and EOS = 0.44), MAE (mean value for SOS = 0.41, LOS = 0.22, and EOS = 0.33), and residual deviance (mean value for SOS = 0.17, LOS = 0.09, and EOS = 0.19) metrics by cross-validation analysis. Overall results suggested that the parameter used for the DFNN model was stable enough to reflect the high accuracy of the model for drought monitoring.
. Figure 5. Cross-validation metrics of deep learning model showing the validation performance for each phenology stage.

Machine Learning Model
Apart from deep learning, two machine learning models, DRF and GBM, were used. These two machine learning algorithms are flexible and robust for classification and regression tasks [16]. DRF is a tree-based algorithm that produces a forest of trees rather than a single tree that uses a randomized subset of candidate features by applying bootstrap sampling to make a final prediction [51]. GBM is also tree-based, using forwardlearning ensemble methods that build a regression tree on all features of the data set for building a predictive model. It works by optimizing the loss function using a weak learner to make a prediction, and by applying an additive model the loss function is minimized. We used the H2O package for DRF and GBM analysis in the R software. To obtain optimum performance, we kept the same value for some important parameters such as 'number of trees', 'max. depth', 'learn_rate', 'sample_rate', fold_assignment = "Random" etc.

Mann-Kendall Test for Drought Trend Analysis
The rank-based Mann-Kendall method, a nonparametric test, was applied to track the increasing and decreasing trends of drought over South Asia at the pixel level during 2001-2016. The Mann-Kendall test, expressed by Mann [76] and explained by Kendall [77], is commonly used for trend detection in the case of time series of environmental, climate and hydrological data [78,79], as well as droughts and aridities of very different spatiotemporal ranges in other very remote regions [80,81]. In general, the Mann-Kendal test follows a monotonic trend rather than a strictly linear trend [82]. The magnitude of the trend was detected by Shen's slope estimator [83] to calculate the upward and downward change of drought levels. A slope with a positive value indicates an increasing trend, while negative values represent a decreasing trend. The equation of Shen's slope is as follows: In Equation (2), β is the trend of SMDI, i and j indicate the interval of the time series, and x i and x j indicate the SMDI value for the year between i and j. The Mann-Kendall test was performed using 'spatialEco' package in R software [84].

Deep Learning and Machine Learning Model Performance with Observed Data and Its Statistical Comparison Using Test Data Set in Detecting Drought Pattern
The discrimination between the deep learning and machine learning approaches for drought distribution in terms of spatial pattern over South Asia is presented in Figure 6. We also observed and compared the spatiotemporal change of drought conditions at three phenology stages during the 2012 drought year, taking drought severity into account. Both deep learning and two machine learning models show relatively low SMDI levels during the SOS and LOS compared to the EOS stage. The northern and southwestern parts of India, Bangladesh, and Nepal had more dry soil during SOS and LOS stages. However, EOS representing drier soil in the southwestern part of India and the northwestern part of Bangladesh. When compared to the observed SMDI distribution, the DFNN and GBM models captured a nearly identical drought pattern. Simultaneously, SMDI levels were underestimated by the DRF model compared to the observed SMDI pattern over the South Asia region. It is evident that at three phenology stages DFNN and GBM caught similar SMDI patterns and well-matched with observed SMDI in terms of spatial distribution over different parts of the South Asia region. To compare further, we statistically evaluated the predicted performance on the test between the deep learning and machine learning models in terms of RMSE, MAE, and MSE (mean square error) during three phenology stages (Table 5). Results presented in Table 5 indicate that DFNN showed better performance with low RMSE, MAE, and MSE than GBM and RF.  As we can see in Figure 7, the frequency distribution of SMDI is almost identical to that of the DFNN and GBM models. Still, less similarity was noticed in frequency distribution for the DRF model. The frequency distribution shows its peak SMDI value of −1.2 and −1.5 for both SOS and LOS stages, respectively. The frequency distribution curve was As we can see in Figure 7, the frequency distribution of SMDI is almost identical to that of the DFNN and GBM models. Still, less similarity was noticed in frequency distribution for the DRF model. The frequency distribution shows its peak SMDI value of −1.2 and −1.5 for both SOS and LOS stages, respectively. The frequency distribution curve was found in the EOS stage with a peak greater than 0 (0.3) for all models. The performance of deep learning and two machine learning models using a tailor diagram (DFNN, DRF, and RF) is presented in Figure 8. The Taylor diagram presented in Figure 8 reveals that DFNN and GBM performed well, demonstrating closer results to each other than the DRF model. Further, the DFNN and GBM model comparison shows that the DFNN performance is slightly better than the GBM model, even though the DFNN and GBM models capture almost similar spatial patterns to the observed data in estimating SMDI. This difference might be due to the multi-layer approach with high tuning parameters that makes computations more practical in the DFNN model than the shallow-depth GBM machine learning model, which has a comparatively low tuning parameter. deep learning and two machine learning models using a tailor diagram (DFNN, DRF, and RF) is presented in Figure 8. The Taylor diagram presented in Figure 8 reveals that DFNN and GBM performed well, demonstrating closer results to each other than the DRF model. Further, the DFNN and GBM model comparison shows that the DFNN performance is slightly better than the GBM model, even though the DFNN and GBM models capture almost similar spatial patterns to the observed data in estimating SMDI. This difference might be due to the multi-layer approach with high tuning parameters that makes computations more practical in the DFNN model than the shallow-depth GBM machine learning model, which has a comparatively low tuning parameter.

Spatiotemporal Pattern of Drought Based on SMDI Using a Deep Learning Model during Three Crop Phenology Stages
To study the effect of precipitation, soil, and vegetation factors on SMDI, we did a spatiotemporal analysis using a deep learning model. Spatially distributed SMDI with the deep learning and two machine learning models using a tailor diagram (DFNN, DRF, and RF) is presented in Figure 8. The Taylor diagram presented in Figure 8 reveals that DFNN and GBM performed well, demonstrating closer results to each other than the DRF model. Further, the DFNN and GBM model comparison shows that the DFNN performance is slightly better than the GBM model, even though the DFNN and GBM models capture almost similar spatial patterns to the observed data in estimating SMDI. This difference might be due to the multi-layer approach with high tuning parameters that makes computations more practical in the DFNN model than the shallow-depth GBM machine learning model, which has a comparatively low tuning parameter.

Spatiotemporal Pattern of Drought Based on SMDI Using a Deep Learning Model during Three Crop Phenology Stages
To study the effect of precipitation, soil, and vegetation factors on SMDI, we did a spatiotemporal analysis using a deep learning model. Spatially distributed SMDI with the

Spatiotemporal Pattern of Drought Based on SMDI Using a Deep Learning Model during Three Crop Phenology Stages
To study the effect of precipitation, soil, and vegetation factors on SMDI, we did a spatiotemporal analysis using a deep learning model. Spatially distributed SMDI with the impact of input model parameters was investigated during three phenology stages: SOS, LOS, and EOS. Figures 9-11 (Figure 12b) demonstrated that drought spatially increases in the central part of Afghanistan and Northwestern part of Pakistan during SOS and LOS stages. Alternatively, Afghanistan, central Pakistan, and the northwestern part of India demonstrated an increasing drought trend during the EOS stage. The erratic rainfall pattern with high potential evapotranspiration causes a significant increasing trend of drought for the region mentioned above. However, the southeastern and southwestern parts of Nepal and northern parts of Bangladesh show a significant decreasing agricultural drought trend due to high rainfall during the pre-monsoon and monsoon seasons. The trend of drought over South Asia during three phenology stages was estimated based on the Mann-Kendal test from 2001-2016, as shown in Figure 12. The spatiotemporal pattern of increasing and decreasing drought trend (Figure 12a) and the corresponding significant level (Figure 12b) demonstrated that drought spatially increases in the central part of Afghanistan and Northwestern part of Pakistan during SOS and LOS stages. Alternatively, Afghanistan, central Pakistan, and the northwestern part of India demonstrated an increasing drought trend during the EOS stage. The erratic rainfall pattern with high potential evapotranspiration causes a significant increasing trend of drought for the region mentioned above. However, the southeastern and southwestern parts of Nepal and northern parts of Bangladesh show a significant decreasing agricultural drought trend due to high rainfall during the pre-monsoon and monsoon seasons. For a clear understanding, we have plotted the observed and predicted SMDI during three phenology stages from 2001-2016 ( Figure 13). The performance of SMDI against the observed SMDI data demonstrates high accuracy prediction. The plots of SMDI are a close to one to one (1:1) line that agrees well with the observed SMDI values. The predicted SMDI looks to be functioning better at the LOS (R2 = 0.52 to 0.94) than the SOS (R2 = 0.57 to 0.90) and EOS (R2 = 0.49 to 0.82) stages. Further, we investigated how much each input variable stimulates the SMDI distribution over South Asia by exploring the variable importance feature to the input function of the deep learning model. Results presented in Table 6 suggested that precipitation factors were the most important factors, contributing significantly to SMDI during the three phenology stages. In addition, among the precipitation factors SPI-6 had the highest contribution effects on SMDI variability simulated by the DFNN model. Further, we investigated how much each input variable stimulates the SMDI distribution over South Asia by exploring the variable importance feature to the input function of the deep learning model. Results presented in Table 6 suggested that precipitation factors were the most important factors, contributing significantly to SMDI during the three phenology stages. In addition, among the precipitation factors SPI-6 had the highest contribution effects on SMDI variability simulated by the DFNN model. In addition, we measured the marginal impact of SPI-6 in three phenology stages presented in Figure 14 with the H 2 O partial plot function of the R program to understand how SPI-6 influences the spatial distribution of SMDI over South Asia. The spatial variability of SMDI followed the variability of SPI-6 to a great extent, indicating that with an increase in SPI-6 value, the SMDI value also increased. The relationship between SPI-6 and SMDI shows a linear trend that describes a unit change of SPI -6 value changes in soil moisture occurs. From the above findings, we can conclude that soil moisture-induced drought increases with a low precipitation rate.  In addition, we measured the marginal impact of SPI-6 in three phenology stages presented in Figure 14 with the H2O partial plot function of the R program to understand how SPI-6 influences the spatial distribution of SMDI over South Asia. The spatial variability of SMDI followed the variability of SPI-6 to a great extent, indicating that with an increase in SPI-6 value, the SMDI value also increased. The relationship between SPI-6 and SMDI shows a linear trend that describes a unit change of SPI -6 value changes in soil moisture occurs. From the above findings, we can conclude that soil moisture-induced drought increases with a low precipitation rate.  We performed a temporal analysis that provided a reasonable estimate of drought severity for each phenology stage that is very sensitive to water stress ( Figure 15) to  33.4%, and 32.2% of the total area, respectively, during the LOS stage. However, comparatively less drought intensity was observed during the EOS stage, influencing 34.64%, 27.48% and 11.2% land of the total area, respectively, in 2002, 2009, and 2014. Furthermore, we also explored the relationship of SMDI simulated by the DFNN model with ground observed SPEI for the evaluation of inter-annual variability of SMDI concerning precipitation. Taking one of the drought-affected stations as an example, the temporal variation of SMDI and SPEI is presented in Figure 16 reveals that SMDI anomaly fluctuations matched well with those of SPEI-3 and SPEI-6 anomaly. Based on the analysis, it can be concluded that SMDI is well correlated with SPEI and SPEI shows a more significant effect on SMDI anomaly. We performed a temporal analysis that provided a reasonable estimate of drought severity for each phenology stage that is very sensitive to water stress ( Figure 15) to understand how much each region is affected by drought. For 33.4%, and 32.2% of the total area, respectively, during the LOS stage. However, comparatively less drought intensity was observed during the EOS stage, influencing 34.64%, 27.48% and 11.2% land of the total area, respectively, in 2002, 2009, and 2014. Furthermore, we also explored the relationship of SMDI simulated by the DFNN model with ground observed SPEI for the evaluation of inter-annual variability of SMDI concerning precipitation. Taking one of the drought-affected stations as an example, the temporal variation of SMDI and SPEI is presented in Figure 16 reveals that SMDI anomaly fluctuations matched well with those of SPEI-3 and SPEI-6 anomaly. Based on the analysis, it can be concluded that SMDI is well correlated with SPEI and SPEI shows a more significant effect on SMDI anomaly.

Uncertainties in SMDI Analysis by DFNN Model
Deep learning-based SMDI prediction is suitable for spatiotemporal drought information that helps government planners to make mitigation and adaptation strategies at the regional level. The deep learning model uses the observation weight column of each incoming data associated with neurons used for bias correction. Thus, the weight of each node determines the output of the entire network. However, the simulated results of drought by the deep learning model mainly depend on learning from past drought information, so an area with no experience of drought is challenging to predict for. Additionally, in our case we used different multi-sensors for remote sensing data with different spatial resolutions, which requires a large number of high-quality training images for better model learning. Moreover, the complexity of the high resolution remote sensing images includes various types of objects with different sizes, rotations, and formats in a single scene, requiring transformation to the same features which creates a problem of robust learning and discriminative illustrations from the objects with deep learning [85]. The above-mentioned causes lead to some uncertainties during the model simulation. For this, we calculated the bias that was generated during the training of the data for each neuron and investigated its spatial pattern ( Figure 17). The bias of all neurons ranges from 0.20-0.55 for three phenology stages (Figure 17a-c). Furthermore, the spatial distribution of bias estimated by calculating the difference between model derived SMDI and observed SMDI presented in Figure 17d-f shows that most of the South Asia region had no bias or tended to zero, but some pixels of the northern and southern parts displayed overestimated and underestimated values compared to the observed values. To minimize bias, we need more training data with a more precise and accurate transformation of remote sensing for better model learning.

Uncertainties in SMDI Analysis by DFNN Model
Deep learning-based SMDI prediction is suitable for spatiotemporal drought information that helps government planners to make mitigation and adaptation strategies at the regional level. The deep learning model uses the observation weight column of each incoming data associated with neurons used for bias correction. Thus, the weight of each node determines the output of the entire network. However, the simulated results of drought by the deep learning model mainly depend on learning from past drought information, so an area with no experience of drought is challenging to predict for. Additionally, in our case we used different multi-sensors for remote sensing data with different spatial resolutions, which requires a large number of high-quality training images for better model learning. Moreover, the complexity of the high resolution remote sensing images includes various types of objects with different sizes, rotations, and formats in a single scene, requiring transformation to the same features which creates a problem of robust learning and discriminative illustrations from the objects with deep learning [85]. The above-mentioned causes lead to some uncertainties during the model simulation. For this, we calculated the bias that was generated during the training of the data for each neuron and investigated its spatial pattern ( Figure 17). The bias of all neurons ranges from 0.20-0.55 for three phenology stages (Figure 17a-c). Furthermore, the spatial distribution of bias estimated by calculating the difference between model derived SMDI and observed SMDI presented in Figure 17d-f shows that most of the South Asia region had no bias or tended to zero, but some pixels of the northern and southern parts displayed overestimated and underestimated values compared to the observed values. To minimize bias, we need more training data with a more precise and accurate transformation of remote sensing for better model learning.

Model Comparison and Performance
This study established the application of a deep learning model for monitoring drought over South Asia. We also used two machine learning models, DRF and GBM, to compare their performance against deep learning. In the case of capturing the spatial pattern of SMDI, the DRF model differs significantly compared to the observed SMDI due to high errors during the training of the data set for each iteration [16,86]. High learning rate and high speed in the case of big data set reduce estimation error for each iteration of the DFNN and GBM models compared to the DRF model, which leads the predicted pattern of SMDI trending towards the observed SMDI pattern. However, the overall performance of deep learning is better than the DRF and RF models because it tends to reduce estimation error with high accuracy due to a multi-layer neural network. Alternatively, DRF requires more computational demand owing to large tree size, which is one of the reasons for low performance over deep learning. Then again, GBM is a slow depth model that offers fewer advantages for dealing with multi-dimensional complex features than the deep learning model [87,88]. These differences might cause the error and low accuracy for predicting SMDI by the DRF and GBM models. Moreover, deep learning has the ability to find the optimal output in the case of high-dimensional data features due to its multi-layer approach. Additionally, many input parameters such as hidden layers with neurons, adaptive learning rate, grid search, dropout function etc. make the computation more practical, more relevant, and advance the underlying algorithms resulting in high predictive accuracy [85]. Sometimes even tree-based non-linear algorithms such as GBM and DRF fail to learn from the data. In that circumstance, deep neural architecture generates non-linear interaction among the variables due to its training stability, facilitating learning the entire network together [75,89].

Drought Variation in Three Phenology Stages
In reality, most of the South Asia region persisted a dry condition and suffered from moderate to severe drought, based on drought classifications presented in Table 1, at three

Model Comparison and Performance
This study established the application of a deep learning model for monitoring drought over South Asia. We also used two machine learning models, DRF and GBM, to compare their performance against deep learning. In the case of capturing the spatial pattern of SMDI, the DRF model differs significantly compared to the observed SMDI due to high errors during the training of the data set for each iteration [16,86]. High learning rate and high speed in the case of big data set reduce estimation error for each iteration of the DFNN and GBM models compared to the DRF model, which leads the predicted pattern of SMDI trending towards the observed SMDI pattern. However, the overall performance of deep learning is better than the DRF and RF models because it tends to reduce estimation error with high accuracy due to a multi-layer neural network. Alternatively, DRF requires more computational demand owing to large tree size, which is one of the reasons for low performance over deep learning. Then again, GBM is a slow depth model that offers fewer advantages for dealing with multidimensional complex features than the deep learning model [87,88]. These differences might cause the error and low accuracy for predicting SMDI by the DRF and GBM models. Moreover, deep learning has the ability to find the optimal output in the case of high-dimensional data features due to its multi-layer approach. Additionally, many input parameters such as hidden layers with neurons, adaptive learning rate, grid search, dropout function etc. make the computation more practical, more relevant, and advance the underlying algorithms resulting in high predictive accuracy [85]. Sometimes even tree-based non-linear algorithms such as GBM and DRF fail to learn from the data. In that circumstance, deep neural architecture generates non-linear interaction among the variables due to its training stability, facilitating learning the entire network together [75,89].

Drought Variation in Three Phenology Stages
In reality, most of the South Asia region persisted a dry condition and suffered from moderate to severe drought, based on drought classifications presented in Table 1, at three phenology stages during the above-mentioned drought years (presented in Section 3.2). Our findings are consistent with the results reported by Mujumdar et al. [90], Neena et al. [91], and Krishnan et al. [92] in which moderate to severe drought in the Indian sub-continent for the years 2002, 2004, 2009, 2014, and 2016 was noticed. The persistence of low soil moisture during the different phenology stages is associated with low precipitation, which brings changes in hydrology and water availability in soil particles. The low precipitation involved with the EI Nino condition over the northwest and northcentral pacific leads to monsoon rainfall deficiency in the South Asia region. Drought most affected the SOS and LOS stages, compared to the EOS stage. SOS is mainly the onset of the monsoon season characterized by a hot summer, which suffers from insufficient rainfall, whereas the LOS stage lasts for the monsoon season with hot and humid conditions, which is also influenced by low rainfall due to Indian Ocean Oscillation and the sea surface temperature gradient between the Indian Ocean and the Bay of Bengal [93,94]. As a result, seasonal drought occurs in SOS and LOS stages. In contrast, the EOS stage, which is the post-monsoon season period, experienced comparatively fewer drought phenomena because the soil particles hold the water that comes from rainfall during the monsoon season. Further, drought gradually expanded from SOS to the LOS stage because crop water requirements reached high levels during LOS due to a rapid increase in the crops' vegetative growth. In addition, high evapotranspiration with rising temperature makes the LOS stage more vulnerable to drought conditions. However, the regional variation of drought trend over South Asia is mainly attributed to the high variability of temperature and rainfall associated with faster changing of the warming phase in the Indian Ocean [95].

SMDI Prediction and Drought Severity Assessment
The SMDI drought index is very robust to short-term dry conditions and very comparable across the season that estimated well the dry conditions in the three phenology stages (SOS, LOS, and EOS), which represents agricultural drought. It is also spatially comparable and varies from one phenology stage to other stages because SMDI is derived from historical value, and the calculation is independent at each grid cell. Though drought affects a large area of each phenology stage, dry and wetness conditions are primarily region-specific, based on precipitation, land cover type, soil texture, and water holding capacity [96][97][98]. The country-level analysis shows the severity of local droughts in India, Pakistan, and Afghanistan during the three phenology stages. This result might be due to the moderate El Niño condition with irregular atmospheric convective activity over the north-central Pacific, which caused rainfall insufficiency over the Indian subcontinent [99][100][101]. Moreover, extended monsoon breaks in the Indian Ocean triggered by ocean-atmosphere dynamical coupling on the intra-seasonal time scale generates drought over the Indian sub-continent [92,102]. The continental-and country-scale analysis aids in understanding the magnitude of local droughts, which may aid in policy formation for post-disaster management to mitigate drought's impact on natural resources.

Conclusions
The present study explored the applicability of a deep learning model to monitor agricultural drought using remote-sensing data. The skill of the model was evaluated using cross-validation data during the training phase of the model. Additionally, the interannual variability of the SMDI was investigated using ground observation data measured with SPEI. The results suggested that the DFNN model was the most-capable tool for monitoring drought across South Asia during three phenology stages. The DFNN model outperformed the other two DRF and GBM models considering precipitation, soil, and vegetation factors. The simulated SMDI by the DFNN model had good consistency with the observed SMDI, and it also matched well with SPEI. The DFNN model estimated yearly drought intensity with high spatial variability across the three phenology stages. The spatial variability was attributed mainly due to the precipitation variability that was examined by the model using relative importance features. Results suggested that droughts in the South Asia region are clearly evident for the years of 2002,2004,2006,2012 27.48% and 11.2% of the total area. The drought prediction system using a deep learning approach facilitates improving agricultural drought monitoring capability over South Asia during crop phenology stages. Moreover, that this study provides a guide of necessity to adopt adaptation strategies or improved management practices by understanding the risk of drought in each crop growing month is one of our contributions. However, the effectiveness and operational ability of the model depend on model input parameters and the size of the training data set. In this regard, incorporation of more drought factors as input parameters with remote sensing-based approaches to meteorological and hydrological drought using different hybrid models could be a future research approach.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data presented in this study can be available on request from the corresponding author.