Temporal and Spatial Comparison of Agricultural Drought Indices from Moderate Resolution Satellite Soil Moisture Data over Northwest Spain

Abstract

During the last decade, a variety of agricultural drought indices have been developed using soil moisture (SM), or any of its surrogates, as the primary drought indicator. In this study, a comprehensive study of four innovative SM-based indices, the Soil Water Deficit Index (SWDI), the Soil Moisture Agricultural Drought Index (SMADI), the Soil Moisture Deficit Index (SMDI) and the Soil Wetness Deficit Index (SWetDI), is conducted over a large semi-arid crop region in northwest Spain. The indices were computed on a weekly basis from June 2010 to December 2016 using 1-km satellite SM estimations from Soil Moisture and Ocean Salinity (SMOS) and/or Moderate Resolution Imaging Spectroradiometer (MODIS) data. The temporal dynamics of the indices were compared to two well-known agricultural drought indices, the atmospheric water deficit (AWD) and the crop moisture index (CMI), to analyze the levels of similarity, correlation, seasonality and number of weeks with drought. In addition, the spatial distribution and intensities of the indices were assessed under dry and wet SM conditions at the beginning of the growing season. The results showed that the SWDI and SMADI were the appropriate indices for developing an efficient drought monitoring system, with higher significant correlation coefficients (R ≈ 0.5–0.8) when comparing with the AWD and CMI, whereas lower values (R ≤ 0.3) were obtained for the SMDI and SWetDI.

1. Introduction

Drought is a complex natural hazard that can affect all climates and produce devastating environmental, economic and social impacts, such as damage to wildlife habitat, diminishment of crop growth or mass migration. While these drought effects have been well documented [1,2], there is not a single and precise definition of drought. It varies among climatic zones and is context-dependent. For this reason, drought assessments rely on the duration, severity and impact area. Droughts are usually classified into four types according to meteorological, agricultural, hydrological and socio-economic criteria [3,4,5,6,7,8,9,10].

Among a wide variety of drought indices, most are climate-based, such as the Standardized Precipitation Index (SPI). The SPI is a meteorological index that compares the precipitation with its multiyear average [11]. The Palmer Drought Severity Index (PDSI) is also a commonly used meteorological index, empirically derived from temperature and precipitation data. The PDSI evaluates the prolonged abnormally wet or dry periods using a two-layer soil model [12]. The PDSI is an effective long-term index, but it does not consider different climate types [13]. To improve this issue, a self-calibrating PDSI (sc-PDSI) was proposed [14]. The Atmospheric Water Deficit (AWD) is also a meteorological index that uses the simple water balance between evapotranspiration (ET) and precipitation [15,16].

To avoid the PDSI limitations [17], the Crop Moisture Index (CMI) was developed. It is an agricultural index that is computed from ET deficits and moisture excess, based on the soil water holding capacity and a subset of calculations required to obtain the PDSI [18]. As an advantage, the CMI responds faster than the PDSI to changing conditions [19]. The two-layer soil model used by Palmer for the PDSI also computes, as an intermediate parameter, the moisture anomaly index (Z-index). The Palmer Z-index is a short-term agricultural index that measures the surface moisture anomaly of the current month without considering antecedent conditions [12,20].

Similar to the PDSI, the Palmer Hydrological Drought Index (PHDI) has been developed. The PHDI uses the same two-layer soil model as the PDSI, but it applies a stricter criterion for determining the ends of drought or wet spells, which results in a longer-term index [12]. In addition, the Soil Water Supply Index (SWSI) was developed as a more appropriate indicator of water availability than the PDSI, especially in regions where rainfall is not the main source of water. Thus, the SWSI is a hydrological index based on the probability distributions of snow runoff, precipitation, reservoir supply, and streamflow on the total surface water resources of the basin [21].

In the literature, there are several drought indicators, i.e., variables or parameters used to describe drought conditions, such as precipitation, temperature, ET, soil moisture (SM) and streamflow [22]. It can be assumed that SM is the governing variable for agricultural drought since this type of drought is considered to begin when the SM availability to plants drops to such a level that it adversely affects the crop yield [7,9]. In addition, SM constrains soil evaporation and plant transpiration in several regions of the world, and it is a storage component for precipitation and radiation anomalies [23]. The soil profile acts as a filter between incoming water and processes that remove water from the hydrological system [24]. Therefore, the SM of the root-zone reflects antecedent conditions and recent precipitation and indicates potential and available water storage [25]. A variety of SM-based agricultural drought indices have been developed over the last decade using SM simulations from models, in situ SM measurements or remotely sensed SM observations.

Some agricultural drought indices have been derived from the probability distributions fitted to a modeled SM. The Variable Infiltration Capacity (VIC) drought index was obtained from the modeled SM through the VIC land surface model [26]. Later, the Soil Moisture Deficit Index (SMDI) was developed from SM anomalies. In this case, the SM was obtained from the Soil Water Assessment Tool (SWAT) hydrological model [27]. The Normalized Soil Moisture (NSM) index was proposed using the SM simulated by the Tiled European Centre for Medium-Range Weather Forecasts (ECMWF) Scheme for Surface Exchanges over Land (TESSEL) model [28]. More recently, the Standardized Soil moisture Index (SSI) and the Multivariate Standardized Drought Index (MDSI) were developed from the modeled SM based on a one-layer water-budget model [29].

All of the previous agricultural drought indices were estimated from hydrological modeling instead of direct root-zone SM measurements. In this regard, the Water Deficit Index (WDI) was introduced as an appropriate way to express SM in terms relative to available water content [30]. Later, the Soil Moisture Index (SMI) was developed from in situ SM data provided by the Automated Weather Data Network (AWDN) and estimated field capacity (FC) and wilting point (WP) water contents [31,32]. Similarly, the Soil Water Deficit Index (SWDI), which is based on the WDI, was proposed using in situ SM data provided by the Soil Moisture Measurements Stations Network of the University of Salamanca (REMEDHUS), and the FC and WP were calculated from water retention curves [33].

During the last decade, surface SM estimations from remote sensing observations have experienced great progress [34,35,36]. Currently, there are two L-band missions in orbit that are specifically devoted to measuring global surface SM [37,38,39,40]. The first mission was launched in November 2009 and is the Soil Moisture and Ocean Salinity (SMOS) satellite from the European Space Agency (ESA) [41]. The second mission was launched in January 2015 and is the Soil Moisture Active Passive (SMAP) satellite from the National Aeronautics and Space Administration (NASA) [42]. In addition, the ESA’s Climate Change Initiative (CCI) generated an SM database using data from passive and active microwave satellite sensors from 1978 to 2015 [43]. These advances allow for the use of satellite surface SM observations for long-term drought monitoring [44].

Preliminary research was performed to analyze the feasibility of using surface SMOS SM data to determine anomalous SM conditions [45]. In this line, the SWDI was applied to SMOS SM data using different methods to estimate FC and WP water contents [46]. The SMOS SM data have also been employed to derive the Soil Moisture Agricultural Drought Index (SMADI) together with the land surface temperature (LST) and the Normalized Difference Vegetation Index (NDVI) from the Moderate Resolution Imaging Spectroradiometer (MODIS) sensor [47,48]. The historical CCI SM data were used to compute the Empirical Standardized Soil Moisture Index (ESSMI) [49]. Alternatively, the Soil Wetness Deficit Index (SWetDI) was developed from the Soil Wetness Index (SWetI) anomaly [50]. In this case, the SWetI played the role of the SM data instead of surface SM observations. The SWetI can be obtained from MODIS LST and NDVI using the LST−NDVI triangle concept to estimate the dry and wet edges [51] following previous interpretations of the LST-NDVI space [52,53]. The quantitative determination of these edges is key for an accurate estimation of ET or SM [54,55]. In addition, the Centre Aval de Traitement des Données SMOS (CATDS) provides the SMOS Drought Index (SDI) at 25 km spatial resolution by including SMOS and MODIS observations into a hydrological model at the root-zone [56].

In drought studies, and particularly in agricultural drought monitoring, the spatial resolution of the data is paramount for the efficient management of drought risks and impacts and water irrigation planning. However, both SMOS and SMAP SM radiometer observations have coarse spatial resolutions (~40 km) that are suitable for global-scale applications [34] but are limited for local and regional drought assessments. During the last several years, there has been a growing interest in developing techniques to enhance the spatial resolution of SMOS and SMAP SM observations. For SMOS, pixel disaggregation approaches are usually based on the synergy of passive microwave and optical visible/infrared (VIS/IR) observations [57,58,59,60,61,62]. Currently, the SMOS Barcelona Expert Centre (BEC) Level 4 (L4) SM maps at 1 km spatial resolution over the Iberian Peninsula are freely distributed (http://bec.icm.csic.es). The downscaling algorithm is based in the LST-NDVI space and integrates SMOS, MODIS and Era Retrospective Analysis (ERA)-Interim data [59,63,64]. An improved downscaling algorithm to obtain those SM maps at 1 km with no size or location limitations is being developed [65].

The goal of this study was to obtain a comprehensive understanding of the spatiotemporal behavior of four innovative agricultural drought indices. To do so, the temporal and spatial performances of the indices were evaluated. This is of special interest in rainfed agro-ecosystems worldwide where the shortages of available water for plant growth directly affect agricultural activities, and these shortages can even induce social and economic conflicts, especially in developing countries. The SM-based indices that were analyzed are SWDI, SMADI, SMDI and SWetDI. All of them use SM (or a surrogate of it) as the main drought indicator. These two facts highlight the novelty of this study: first, the use of the SM variable, which is relatively new in drought studies; and, second, the agricultural drought perspective instead of the more common meteorological approach. The indices were computed over northwest Spain in the Castilla y León region on a weekly basis from June 2010 to December 2016 using the SMOS BEC L4 SM and/or MODIS reflectance and LST data at 1 km spatial resolution. In the case of SWDI, two different approaches to estimate the FC and WP were separately assessed. The temporal dynamics of the indices were compared to two well-known agricultural drought indices: the AWD and the CMI. The spatial distribution of each index and its intensity were assessed by comparing the estimated drought index maps under dry and wet SM conditions at the beginning of the growing season. This research could help enlighten the choice of the more appropriate drought index for developing an accurate agricultural drought monitoring system.

2. Data and Methods

2.1. In Situ Data

Meteorological data from the Soil Moisture Measurement Stations Network of the University of Salamanca (REMEDHUS) [66]; Inforiego (www.inforiego.org), which is managed by the Agriculture Technological Institute of Castilla y León (ITACYL); and the Spanish Meteorological Agency (AEMet, www.aemet.es) weather station networks were used in this work. These three networks are in the Castilla y León region (northwest Spain), which is one of the largest agricultural regions in the European Union (Figure 1, top right). The study area is a large plain that is approximately 65,000 km² with an average altitude of 800 m above sea level and is surrounded by mountains to the north, south and east. The area has a continental semi-arid Mediterranean climate with a mean annual temperature of 11.8 °C and mean annual precipitation of 450 mm [67]. The principal land use is agriculture, including cereals (73.5%), industrial crops (13%), forages (8%), legumes (4.5%), tubers (0.7%) and vegetables (0.3%) [68]. A land uses map at a spatial resolution of 25 m from the ITACYL map service (http://mcsncyl.itacyl.es/en/descarga) is displayed in Figure 1, where the different crops and natural areas were grouped into seven categories (agricultural, artificial, rocky/bare, pasture, scrubland, forest and water bodies).

In this study, the meteorological variables (i.e., air temperature, relative humidity, precipitation, solar radiation and wind speed) measured by a REMEDHUS weather station (Villamor) and 22 Inforiego agro-meteorological stations were used from June 2010 to December 2016. In the case of AEMet, long-term temperature and precipitation data provided by 6 stations from 1985 to 2016 were used. In addition, a soil database of the Duero River basin with information on sand, clay, silt and organic matter (MO) contents, provided by the ITACYL (http://suelos.itacyl.es/) was used to obtain soil texture and organic content maps and to derive specific pedotransfer functions (PTFs) for this region. The ITACYL database consists of approximately ten thousand punctual soil samples.

2.2. Satellite Data

The SMOS BEC L4 SM v.3 maps from June 2010 to December 2016 at 1 km spatial resolution were used. They were obtained using a downscaling algorithm, which was applied separately for both daily SMOS ascending (6:00 a.m.) and descending (6:00 p.m.) orbits. This algorithm combines daily SMOS brightness temperature (Level 1C v.620) and SM (Level 2 v.620) data, 16-day Terra MODIS NDVI (MOD13A2) data and daily ERA-Interim LST data [59,63,64].

The most recent daily Aqua MODIS surface reflectance (SR) data in bands B1 (~660 nm) and B2 (~860 nm) at 500 m spatial resolution (MYD09GA v.6) and MODIS LST data at 1 km spatial resolution (MYD11A1 v.6) from June 2010 to December 2016 were employed in this study. Aqua was used instead of Terra because the Aqua day-time LST (at 1:30 p.m. of acquisition local time) is closer to the daily maximum LST and was more correlated to the SM than the Terra day-time LST (at 10:30 a.m.) [63,69]. The SR is already corrected for atmospheric conditions such as gases, aerosols and Rayleigh scattering. The nominal accuracy of the LST is 1 °C under clear sky conditions. Both products are projected to a tile-based sinusoidal grid. Four tiles covering the study area (h17v04, h17v05, h18v04 and h18v05) were selected.

2.3. Data Processing

To analyze the drought events and patterns of different agricultural drought indices over the Castilla y León region, the agricultural areas were delimited. In this region, crop areas generally occur in the central part, out of the mountains. Thus, the first step consisted of applying an agricultural mask to the entire region based on the digital elevation model (DEM) from the ITACYL map service (Figure 1). This DEM was resampled from the 25-m grid to a regular 1-km grid in the World Geodetic System (WGS)-84 datum. The agricultural mask preserves areas with altitudes lower than 1100 m. Hence, drought indices were only estimated at this spatial domain, and all daily input data needed to compute the indices were previously filtered.

In addition, the soil data were also filtered to discard samples over forested regions, and the soil samples with MO contents higher than 4% were removed. The remaining 9951 samples were interpolated using the inverse distance weighting to the regular 1-km grid in WGS-84 and then filtered with the agricultural mask. The resulting maps of sand, clay, silt and MO contents at 1 km over the agricultural areas were used (Figure 2).

The SMOS BEC L4 SM v.3 maps were first clipped to select a rectangular area of the Castilla y León region (39.85°N–43.35°N; 7.35°W–1.65°W) and resampled to the regular 1-km grid in WGS-84. Later, the ascending and descending maps were averaged to obtain a unique map for each day, and then the daily maps were filtered with the agricultural mask. Finally, these daily maps were weekly averaged. The pixels with less than four SM measurements within the week were discarded.

Regarding the MODIS products, for each day, the four tiles of MODIS SR and LST were mosaicked and clipped to select the rectangular area over the Castilla y León region. These maps were resampled from the sinusoidal grid to the regular grid in WGS84 at their original spatial resolutions (SR at 500 m and LST at 1 km) using the automated batch processing of the MODIS Reprojection Tool (MRT) [70]. The SR in bands B1 and B2 was filtered using their corresponding state (internal cloud algorithm, MOD35 snow/ice and internal snow mask) and quality assurance flags (highest quality at B1 and B2). Then, the daily MODIS NDVI at 500 m was computed as follows:

$$\mathrm{NDVI}=\frac{{\mathrm{SR}}_{\mathrm{B}2}-{\mathrm{SR}}_{\mathrm{B}1}}{{\mathrm{SR}}_{\mathrm{B}2}+{\mathrm{SR}}_{\mathrm{B}1}}$$

(1)

where SR_B1 and SR_B2 correspond to the surface reflectance at bands B1 and B2, respectively. NDVI pixels at 500 m were upscaled to the same regular 1-km grid of the LST using a simple average. Both daily LST and NDVI at 1 km were filtered with the agricultural mask. Finally, the LST and NDVI maps were weekly averaged.

2.4. Estimation of Drought Indices

2.4.1. Atmospheric Water Deficit (AWD)

The daily ET was estimated by the Penman–Monteith method [71], and the measured precipitation from the Inforiego and Villamor stations (see locations in Figure 1) were used. The AWD was proposed as a 7-day running sum of ET minus a 7-day running sum of precipitation [16]. In this study, the weekly AWD was reversely calculated as:

$${\mathrm{AWD}}_i={\mathrm{P}}_i-{\mathrm{ET}}_i$$

(2)

where P_i corresponds to the sum of precipitation and ET_i denotes the sum of evapotranspiration for the current week i.

2.4.2. Crop Moisture Index (CMI)

The daily air temperature and precipitation from 6 AEMet stations over the Castilla y León region (see locations in Figure 1) from 1985 to 2016 were employed to compute the weekly CMI using the sc-PDSI tool [72]. The particular inputs required for each station were: (i) the mean monthly (and weekly) air temperature for each of the 12 months (52 weeks) of the year; (ii) the accumulated monthly (and weekly) precipitation for each of the 12 months (52 weeks) of the year; (iii) the mean monthly (and weekly) air temperature over all the years; and (iv) the total available water (TAW).

The TAW in the root-zone was computed with the following equation [71]:

$$\mathrm{TAW}=1000\times \left(\mathrm{FC}-\mathrm{WP}\right)\times {\mathrm{z}}_{\mathrm{r}}$$

(3)

where FC and WP are the soil water content at field capacity and wilting point, respectively, and z_r stands for the effective rooting depth. To estimate the TAW, an effective root length of 1 m was taken into account, in agreement with the maximum effective rooting depth for the most common cereals in the area (wheat, barley and oats) [71].

The soil water parameters (FC and WP) were estimated from the sand, clay and MO maps at 1 km (see maps in Figure 2) using the value of the pixel covering each AEMet station and the PTFs proposed by Rawls et al. [73]:

$${\mathrm{FC}}_{\mathrm{Rawls}}=0.2576-0.0020\times \mathrm{S}+0.0036\times \mathrm{C}+0.0299\times \mathrm{MO}$$

(4)

$${\mathrm{WP}}_{\mathrm{Rawls}}=0.0260+0.0050\times \mathrm{C}+0.0158\times \mathrm{MO}$$

(5)

where S is sand, C is clay, and MO corresponds to the organic matter contents (all in percentages).

2.4.3. Soil Water Deficit Index (SWDI)

Considering that only a fraction of the soil water content is available for plant growth (available water content, AWC), the SWDI characterizes the agricultural drought directly from the SM referred to as the AWC, which is the difference between the FC and WP. The weekly SWDI was computed as [33,46]:

$${\mathrm{SWDI}}_i=\left(\frac{{\mathrm{SM}}_i-\mathrm{FC}}{\mathrm{AWC}}\right)\times 10$$

(6)

where SM_i is the SMOS BEC L4 SM at 1 km for the current week i.

The SWDI was formulated to obtain a range of values with agricultural meaning. Its interpretation is related to the p factor of readily available water (RAW), which is the average fraction of TAW that can be depleted from the root-zone before moisture stress occurs [71].

To analyze the impact of different PTFs on the SWDI estimation, two different expressions for the same soil water parameters were used. Firstly, the FC and WP were estimated using Equations (4) and (5), respectively. Alternatively, a multiple regression analysis was performed using the soil water parameters as the dependent variables and the soil database (sand, clay, silt and organic matter contents) as the independent variables to obtain specific PTFs for the Castilla y León region. The best PTFs that were obtained, in terms of coefficient of determination and bias, were as follows:

$${\mathrm{FC}}_{\mathrm{CYL}}=0.41-0.004\times \mathrm{S}+0.01\times \mathrm{MO}$$

(7)

$${\mathrm{WP}}_{\mathrm{CYL}}=0.002+0.005\times \mathrm{C}+0.015\times \mathrm{MO}$$

(8)

Therefore, two possible SWDI approaches were estimated: (i) the SWDI_Rawls, which employed the PTFs of Rawls et al. in Equations (4) and (5); and (ii) the SWDI_CYL, which used the specific PTFs for the Castilla and León region in Equations (7) and (8).

2.4.4. Soil Moisture Agricultural Drought Index (SMADI)

Three condition indices are needed to compute the SMADI: (i) the Soil Moisture Condition Index (SMCI), (ii) the Modified Temperature Condition Index (MTCI), and (iii) the Vegetation Condition Index (VCI). They were calculated for each week using the following expressions [47]:

$${\mathrm{SMCI}}_i=\frac{{\mathrm{SM}}_{\max }-{\mathrm{SM}}_i}{{\mathrm{SM}}_{\max }-{\mathrm{SM}}_{\min }}$$

(9)

$${\mathrm{MTCI}}_i=\frac{{\mathrm{LST}}_i-{\mathrm{LST}}_{\min }}{{\mathrm{LST}}_{\max }-{\mathrm{LST}}_{\min }}$$

(10)

$${\mathrm{VCI}}_i=\frac{{\mathrm{NDVI}}_i-{\mathrm{NDVI}}_{\min }}{{\mathrm{NDVI}}_{\max }-{\mathrm{NDVI}}_{\min }}$$

(11)

where SM_i is the SMOS BEC L4 SM at 1 km, LST_i corresponds to the MODIS LST at 1 km, and NDVI_i denotes the MODIS NDVI at 1 km, all of them for the current week i. The subscripts max and min denote the maximum and minimum values of the pixel for the entire study period, respectively.

Then, the weekly SMADI was computed as [47,48]:

$${\mathrm{SMADI}}_i={\mathrm{SMCI}}_i\times \frac{{\mathrm{MTCI}}_i}{{\mathrm{VCI}}_{i+1}}$$

(12)

where SMCI_i is the SMCI of the current week i, MTCI_i stands for the MTCI of week i, and VCI_i+1 corresponds to the ensuing week i + 1.

2.4.5. Soil Moisture Deficit Index (SMDI)

The first step to calculate the SMDI was to obtain the SMOS soil moisture deficit data (SD_SMOS, in percentage), which were estimated for each week i and year k as a weekly SM anomaly with respect to the median value, using the following conditional expressions [27]:

$${\mathrm{SD}}_{\mathrm{SMOS},k,i}=\frac{{\mathrm{SM}}_{k,i}-{\mathrm{SM}}_{\mathrm{median},i}}{{\mathrm{SM}}_{\mathrm{median},i}-{\mathrm{SM}}_{\min ,i}}\times 100{\mathrm{if}\mathrm{SM}}_{k,i}\le {\mathrm{SM}}_{\mathrm{median},i}$$

(13)

$${\mathrm{SD}}_{\mathrm{SMOS},k,i}=\frac{{\mathrm{SM}}_{k,i}-{\mathrm{SM}}_{\mathrm{median},i}}{{\mathrm{SM}}_{\max ,i}-{\mathrm{SM}}_{\mathrm{median},i}}\times 100{\mathrm{if}\mathrm{SM}}_{k,i}>{\mathrm{SM}}_{\mathrm{median},i}$$

(14)

where SM_k,i is the SMOS BEC L4 SM for the current week i (1 to 52) and year k (2010 to 2016), and the SM_median,i, SM_max,i and SM_min,i stand for the long-term median, maximum and minimum values of the current week i during the entire study period, respectively.

Then, the weekly SMDI was computed as [27]:

$${\mathrm{SMDI}}_i=0.5\times {\mathrm{SMDI}}_{i-1}+\frac{{\mathrm{SD}}_{\mathrm{SMOS},i}}{50}$$

(15)

where SMDI_i−1 is the SMDI of the previous week i − 1 and SD_SMOS,i corresponds to the SD_SMOS of the current week i. When data were not available from the previous week, the nearest antecedent value of SMDI was chosen. The SMDI was initialized as SMDI₁ = SD₁/50.

2.4.6. Soil Wetness Deficit Index (SWetDI)

Instead of using SM observations, the SWetDI used a parameter, the Soil Wetness Index (SWetI), which mimics the role of SM. The weekly SWetI was estimated as [51]:

$${\mathrm{SWetI}}_{i,n}=\frac{{\mathrm{T}}_{\max ,i,n}-{\mathrm{Ts}}_{i,n}}{{\mathrm{T}}_{\max ,i,n}-{\mathrm{T}}_{\min ,i,n}}$$

(16)

where Ts_i,n is the MODIS LST at 1 km for week i and pixel n, and T_max,i,n and T_min,i,n correspond to the maximum and minimum values of observed MODIS LST at 1 km for each n-th MODIS NDVI at 1 km and week i. These values were derived weekly from the daily observations of the LST-NDVI space at the agricultural domain (Figure 3) using linear relations between LST and NDVI on the dry and wet edges:

$${\mathrm{T}}_{\max ,n}=a\times {\mathrm{NDVI}}_n+b,\mathrm{dry}\mathrm{edge}$$

(17)

$${\mathrm{T}}_{\min ,n}=c\times {\mathrm{NDVI}}_n+d,\mathrm{wet}\mathrm{edge}$$

(18)

where a and c are the slopes, and b and d are the offsets of the diagonal lines on the dry and wet edges, respectively. These coefficients were computed for each week using robust linear regressions [74] by applying an iterative algorithm similar to a previous study [54].

The MODIS soil moisture deficit (SD_MODIS, in percentage) has a similar expression as the previous SD_SMOS but uses the SWetDI anomaly. It was estimated for each week i and month j as a monthly SWetI anomaly with respect to the mean value by the following expression [50]:

$${\mathrm{SD}}_{\mathrm{MODIS},j,i}=\frac{{\mathrm{SWetI}}_{j,i}-{\mathrm{SWetI}}_{\mathrm{mean},j}}{{\mathrm{SWetI}}_{\max }-{\mathrm{SWetI}}_{\min }}\times 100$$

(19)

where SWetI_j,i stands for the SWetI of the current week i (1–52) and month j (1–12), SWetI_median,j is the long-term mean SWetI for month j, and SWetI_max and SWetI_min correspond the maximum and minimum SWetI for the entire study period, respectively.

The weekly SWetDI was estimated, similar to the SMDI, as [50]:

$${\mathrm{SWetnessDI}}_i=0.5\times {\mathrm{SWetnessDI}}_{i-1}+\frac{{\mathrm{SD}}_{\mathrm{MODIS},i}}{50}$$

(20)

where SWetnessDI_i−1 is the SWetDI of the previous week i − 1 and SD_MODIS,i corresponds to the SD_MODIS of the current week i. When data were not available from the previous week, the nearest antecedent value of SWetDI was chosen. The SWetDI was also initialized as SWetDI₁ = SD₁/50.

2.5. Comparison Strategy

In this study, weekly comparisons were made because one week is the time interval commonly used by farmers for irrigation schedules [15]. Therefore, weekly maps of the drought indices (SWDI_Rawls, SWDI_CYL, SMADI, SMDI and SWetDI) were estimated over the agricultural areas of the Castilla y León region. The time-series of the pixels overlapping each in situ station were extracted and compared with the AWD time-series computed at the Inforiego and Villamor stations and the CMI time-series computed at the AEMet stations. For these comparisons, Pearson’s correlation coefficient (R) was applied. It is obvious that in situ and satellite observations have different spatial resolution, but this validation strategy has been widely used in literature, where moderate and even low resolution satellite drought indices were compared to punctual data-based drought indices [46,47,48,75,76,77,78,79]. Besides, after an additional test (results not shown) consisting of computing again the correlation between the SMDI and SWetDI with AWD and CMI, after discarding the first six months of the time-series, it was verified that the value of initialization in the SMDI and SWetDI does not affect results and neither conclusions obtained from this study. The number of weeks when drought occurred was also computed for each drought index. The occurrence of drought depends on the established thresholds of each index according to their categories (

Table 1. Categories of each drought index.

		AWD	CMI	SWDI	SMADI	SMDI	SWetDI
Dynamic Range		−∞ to +∞	−∞ to +∞	−∞ to +∞	0 to +∞	−4 to +4	−4 to+ 4
No drought		0 or more	0 or more	0 or more	0 to 1	0 or more	0 or more
Drought	Mild	less than 0	−2 to −0.01	−2 to −0.01	1.01 to 2	−1 to −0.01	−1 to −0.01
	Moderate		−3 to −2.01	−5 to −2.01	2.01 to 3	−2 to −1.01	−2 to −1.01
	Severe		less than −3	−10 to −5.01	3.01 to 4	−3 to −2.01	−3 to −2.01
	Extreme		-	less than −10	more than 4	−4 to −3.01	−4 to −3.01

). In this study, all intensities of drought (mild, moderate, severe and extreme) were considered to calculate the number of drought weeks.

The occurrence of drought of each index was assessed by several categorical statistics referred to an observed event (

Table 2. Categorical statistics, equations, dynamic ranges and perfect scores.

Categorical Statistic	Equation	Dynamic Range	Perfect Score
Probability of Detection (POD)	$\mathrm{POD}=\frac{\mathrm{A}}{\mathrm{A}+\mathrm{C}}$	0 to 1	1
Probability of False Detection (POFD)	$\mathrm{POFD}=\frac{\mathrm{B}}{\mathrm{B}+\mathrm{D}}$	0 to 1	0
False Alarm Ratio (FAR)	$\mathrm{FAR}=\frac{\mathrm{B}}{\mathrm{A}+\mathrm{B}}$	0 to 1	0
Frequency Bias (FB)	$\mathrm{FB}=\frac{\mathrm{A}+\mathrm{B}}{\mathrm{A}+\mathrm{C}}$	−∞ to +∞	1

A = hits; B = false alarms; C = misses; and D = correct negatives

): (i) the probability of detection (POD); (ii) the probability of false detection (POFD); (iii) the false alarm ratio (FAR); and (iv) the frequency bias (FB). These statistics are widely used for weather forecast verification in which the forecast is compared with an observed event or reference at four possible cases based on the A, B, C and D values of a contingency table [80,81]. In this study, the AWD and the CMI were assumed to be the observed events. Thus, A corresponds to the hits, i.e., the weeks when both the reference and the index indicate drought. B denotes the false alarms, which are the weeks when the index indicates drought but the reference indicates no drought. C stands for the misses, which are the weeks when the reference indicates drought but the index indicates no drought. D denotes the correct negatives, i.e., the weeks when both the reference and the index indicate no drought.

According to this contingency table, the POD denotes what fraction of the observed drought weeks was correctly estimated. The POFD computes what fraction of the observed no drought weeks was incorrectly estimated as drought. The FAR measures what fraction of the estimated drought weeks did not occur. The FB expresses the systematic differences between the reference and the index using the frequency of the number of drought weeks. All equations used to calculate these categorical statistics, as well as their dynamic ranges and perfect scores, are included in Table 2 [82]. These statistics were displayed in boxplots, where the central mark indicates the median, the bottom edge of the box displays the 25th percentile, and the top edge shows the 75th percentile of all stations. The whiskers extend to the most extreme data points that were not considered outliers, and the outliers are plotted individually using red crosses.

The last comparison assessed the spatial distribution of the estimated drought maps for each index. Thus, the drought indices were spatially compared at two representative weeks of the study period: (i) 8–14 April 2012; and (ii) 9–15 April 2013. For most crops of the Castilla y León region, the plant growing season lasts from approximately April to October. The two selected weeks correspond to the beginning of the growing season during two years with dry and wet SM conditions, respectively. The meteorological drought was very intense in 2012 [83], and the SM content was very low in this region. Conversely, in 2013, the amount of precipitation was higher than the long-term mean. The second reason to choose these specific weeks relies on minimizing the presence of clouds. Since optical VIS/IR data are affected by clouds, there are some gaps in the daily MODIS LST and NDVI observations. However, as done for the other datasets, the daily MODIS maps were weekly averaged, reducing these data gaps. The indices that could be affected by them were SMADI and SWetDI. Furthermore, the chosen weeks have almost the maximum spatial coverage, with only a few pixels with no data. The SMOS BEC L4 SM data were not affected by this issue because its downscaling algorithm employs a filled composite of 16-day NDVI at 1 km.

The flowchart in Figure 4 graphically displays all the data processing and the methodology applied to compute the analyzed drought indices of this study.

3. Results and Discussion

3.1. Time-Series Comparison

The weekly evolution of the AWD, the CMI and the drought indices at the two Inforiego and AEMet stations is shown in Figure 5 and Figure 6, respectively. VA01 and VA101 were selected as examples from the Inforiego network because they correspond to the wettest and driest stations along the entire study period, with the SM varying from 0.03 to 0.40 m³/m³ and from 0.02 to 0.27 m³/m³, respectively. In the AEMet network, the same criterion was applied to select the Zamora and León stations, with a SM variability from 0.03 to 0.39 m³/m³ and from 0.03 to 0.28 m³/m³, respectively. Although the SMADI was defined in the opposite sense with respect to the other drought indices, for a better comparison, its vertical axis was oriented downward.

It was clearly observed that the SWDI_Rawls, SWDI_CYL and SMADI displayed similar seasonal cycles to the AWD and CMI (Figure 5a–f and Figure 6a–f). In contrast, the SMDI and SWetDI did not exhibit any seasonal patterns (Figure 5g–j and Figure 6g–j). This result can be a consequence of using the SM deficit obtained from the SM or SWetI anomalies instead of directly using the SM or SWetI for the calculation of the two indices, since these anomalies remove the inherent seasonality of the original series of SM or the SWetI. By definition, an agricultural drought index should reflect the available soil water conditions for vegetation growth. This seasonal cycle was effectively captured by the AWD, CMI, SWDI and SMADI, being a distinctive feature of these drought indices.

The temporal evolution of the AWD and both SWDI approaches were nearly synchronized (Figure 5a–d). The one-week lag observed between the AWD and the SWDI computed from the in situ SM data at the root-zone in a previous study [33] was not observed when using the SMOS-derived SWDI [46,78]. This is because the SMOS SM is representative of only the top 5 cm of the soil profile, where soil water evaporation preferably occurs and water moves faster. Moreover, the SMOS SM data have a quicker response to rainfall events and dry-downs than the in situ SM data, even at the same soil layer [66,67].

When comparing the CMI and both SWDI approaches, there were some drought breaks during the dry periods when the CMI was zero, which were not detected by the SWDI or the other indices (Figure 5a–d). This fact was previously observed, and it is related to small rainfall events that do not have effective influences on the recharge of soil water [33,46,47].

The SMADI was delayed respect to the AWD at the onset of the drought events (Figure 5e,f). There is a well-known time-lag effect in the occurrence of a drought and the changes in NDVI [75,84,85,86,87]. Previous studies revealed a delay of 20 days, 16–24 days or even one month between a rainfall event and changes in NDVI [88,89,90]. This may suggest that the one-week delayed NDVI could be not long enough for the SMADI computation. However, a lag of 10 days was found between the NDVI anomalies and soil moisture in grasslands, croplands and shrublands [91], which is in agreement with the one-week delayed used in this work.

As expected, the SWDI_Rawls and SWDI_CYL had the same temporal evolution and only differed in their absolute values and dynamic ranges (Figure 5a–d and Figure 6a–d). These differences were produced by the different PTFs used to estimate the soil water parameters. Both SWDI approaches estimated drought in VA101 and León in almost the entire study period, especially the SWDI_Rawls (Figure 5b,d and Figure 6b,d). This result indicated that the SWDI is very sensitive to the accuracy of the FC and WP, particularly when the SM is very low. In addition, the accuracy of any method based on PTFs depends on the robustness of the soil database, as previously found [46].

The SMADI exhibited the smoothest temporal evolution, i.e., the signals of the other indices were noisier than the signals of the SMADI, but there were few peaks with very high SMADI values at particular weeks. It is highlighted that the SMADI displayed higher values in León than in Zamora since León is drier than Zamora in terms of SM (Figure 6e,f). Nevertheless, VA101 was drier than VA01, and this effect was not observed. In contrast, the SMADI values were higher in VA01 than in VA101, which means that the SMADI estimated a stronger drought in VA01 than in VA101 (Figure 6e,f). This result can be explained by the integration of the LST and the NDVI additional information in SMADI [47].

The SMDI showed the highest temporal variability with sudden events of drought, changing from non-drought to drought and vice versa during very short time intervals, i.e., going upward and downward between the maximum and minimum limits (Figure 5g,h and Figure 6g,h). This result may be due to the SMDI was originally developed from an SM time-series of 98 years [27], which was much larger than the study period (6 and a half years). Since the period is much shorter, the SM deviation of a given week with respect to its weekly median along the study period is expected to be higher. Despite this, and considering the SMDI trend, most weeks of the SMDI displayed drought conditions in 2012. This result agrees with the exceptional meteorological drought of 2012 [83]. The SMDI displayed no drought in most weeks of 2013. No significant differences were observed between the SMDI at the wettest or the driest stations.

The SWetDI had a very low dynamic range, usually from −1 to +1 and did not reach its minimum or maximum limits (−4 to +4) in any week of the study period (Figure 5i,j and Figure 6i,j). This variability of the SWetDI was similar to that obtained in a previous study, where the SWetDI was calculated over Iran, a region with an annual precipitation of 130 mm per year (much drier than the Castilla y León region), and a variability ranging from −1.7 to +1.2 in a rainfed farmland location was displayed [50]. Nevertheless, the size of the observed scene and the different climates or land uses could have an important impact on obtaining the dry and wet edges through the LST-NDVI space, as well as in the estimation of SWetI. Furthermore, the SWetDI trend was not capable of capturing the drought event of 2012. As in the case of the SMDI, no significant differences were observed between the SWetDI at the wettest or the driest stations. Similar results were obtained from the time-series of drought indices at the other Inforiego and Villamor stations, or at the other AEMet stations (not shown).

3.2. Correlation Analysis

For each Inforiego and Villamor station, the correlation coefficients of the AWD and the drought indices are shown in

Table 3. Pearson’s correlation coefficients (R) of weekly Atmospheric Water Deficit (AWD) and Soil Water Deficit Index (SWDI)_Rawls, SWDI_CYL, Soil Moisture Agricultural Drought Index (SMADI), Soil Moisture Deficit Index (SMDI) or Soil Wetness Deficit Index (SWetDI) at the Inforiego and Villamor stations.

Station	SWDIRawls	SWDICYL	SMADI	SMDI	SWetDI
AV01	0.76	0.76	−0.71	0.18	0.02 *
BU03	0.78	0.78	−0.51	0.19	0.09 *
BU04	0.68	0.68	−0.62	0.19	0.15
BU05	0.76	0.76	−0.51	0.18	0.14
LE03	0.69	0.69	−0.41	0.20	0.15
LE04	0.71	0.71	−0.19	0.20	0.06 *
LE08	0.78	0.78	−0.50	0.28	0.14
P02	0.74	0.74	−0.55	0.20	0.11 *
P04	0.77	0.77	−0.60	0.25	0.11 *
P06	0.71	0.71	−0.40	0.24	0.05 *
SA101	0.74	0.74	−0.61	0.20	0.08 *
SA102	0.78	0.78	−0.53	0.19	0.08 *
SG02	0.74	0.74	−0.77	0.16	0.11 *
SO02	0.79	0.79	−0.62	0.17	0.13
VA01	0.80	0.80	−0.72	0.24	0.15
VA02	0.77	0.77	−0.64	0.20	0.09 *
VA05	0.80	0.80	−0.68	0.17	0.11 *
VA06	0.77	0.77	−0.29	0.23	0.09 *
VA08	0.78	0.78	−0.67	0.21	0.07 *
VA101	0.75	0.75	−0.83	0.14	0.00 *
ZA02	0.77	0.77	−0.19	0.19	0.06 *
ZA05	0.75	0.75	−0.67	0.22	0.13
VILLAMOR	0.79	0.79	−0.69	0.26	0.12 *

* Non-significant correlations (ρ_val > 0.05).

. Note that the correlation between the AWD and SWDI_Rawls was equal to the correlation between the AWD and SWDI_CYL at all stations (R ≈ 0.68–0.80). This result indicates that there was a very high agreement between the AWD and both SWDI approaches, regardless of the PTFs used, and the SWDI adequately reproduced the soil water balance dynamics. The PTFs of Rawls et al. [73] worked as well as the specific PTFs for the Castilla y León region when estimating the soil water parameters. Therefore, it can be concluded that the SWDI can be calculated in any region with available soil data by applying a suitable standard and adequately tested PTFs. The correlations were similar to those obtained in previous research [46], where the AWD was compared with a SMOS-derived SWDI (R = 0.83) using estimations of FC and WP from soil samples in a laboratory.

There was a negative correlation between the AWD and SMADI, which showed a large variation among stations (R ≈ −0.19 to −0.83). Some stations (BU03, BU05, LE03, LE04, LE08, P06, SA102, VA06 and ZA02) showed a correlation weaker than −0.55. Although these stations are located over rainfed cereals plots, they are surrounded by areas of irrigated crops. This fine-scale detail was hidden in the scale of the remote sensing observations, which represent an integrated value over an area. Thus, it is likely that the response of the irrigated crops was mixed with the response of the rainfed stations. The correlations at the rest of the stations were R ≈ −0.60 to −0.83, which indicates a high agreement between AWD and SMADI.

The correlation between the AWD and SMDI was very low in all stations (R ≤ 0.28). Similarly, the correlation between the AWD and SWetDI was very low in all stations (R ≤ 0.15) after the non-significant values were discarded (ρ_val > 0.05). Not surprisingly, this could be due to the absence of the seasonal dynamics in both the SMDI and SWetDI time-series, as previously shown in Figure 5g–j. Another possible hypothesis about the poor results of SWetDI could be related to the physical interpretation of the LST-NDVI triangle, because recent research has pointed out better results of the LST-NDVI approach when decoupling the soil evaporation and vegetation transpiration [55].

The correlation coefficients between the CMI and the other SM-based drought indices for each AEMet station are shown in

Table 4. Pearson’s correlation coefficients (R) of weekly Crop Moisture Index (CMI) and Soil Water Deficit Index (SWDI)_Rawls, SWDI_CYL, Soil Moisture Agricultural Drought Index (SMADI), Soil Moisture Deficit Index (SMDI) or Soil Wetness Deficit Index (SWetDI) for the Spanish Meteorological Agency (AEMet) stations.

Station	SWDIRawls	SWDICYL	SMADI	SMDI	SWetDI
BURGOS	0.60	0.60	−0.58	0.44	0.25
LEÓN	0.69	0.69	−0.61	0.39	0.19
SALAMANCA	0.66	0.66	−0.57	0.40	0.08 *
SORIA	0.48	0.48	−0.54	0.41	0.28
VALLADOLID	0.70	0.70	−0.32	0.42	0.07 *
ZAMORA	0.69	0.69	−0.69	0.44	0.13

* Non-significant correlations (ρ_val > 0.05).

. As previously observed with the AWD, the same correlation was obtained between the CMI and both SWDI approaches (R ≈ 0.48–0.70). This result is similar (R = 0.72) to that obtained in a previous study [46]. Nevertheless, the agreement between the SWDI and CMI was lower than the agreement with AWD, as previously found [46].

The correlation between the SMADI and CMI (R ≈ −0.32 to −0.69) was similar to the value obtained in a previous study (R ≈ −0.71) [47]. This correlation between SMADI and CMI was also weaker than the correlation with AWD. This may be because the CMI produced several near-zero values, both at the beginning and at the end of the growing season [47]. The correlation between the SMDI and CMI was weak (R ≈ 0.39–0.44), and the correlation between the CMI and SWetDI was very weak (R ≤ 0.28). As in the case of AWD, this result could be related to the lack of temporal patterns in both the SMDI and SWetDI time-series, as previously shown in Figure 6g–j.

3.3. Drought Weeks Captured

The number of drought weeks captured by each index is shown in Figure 7. Since each station had a different number of weeks with available data along the study period, the comparison was made using the percentage of the total number of weeks. On the one hand, the AWD identified between 75.2% and 88.4% of the drought weeks (Figure 6a). On the other hand, the CMI captured between 31.7% and 45.3% (Figure 6b). These results were similar to those obtained in a previous study in the REMEDHUS area (June 2010–December 2014), where 73.1–84.6% and 32.7–57.7% of the drought weeks per year were estimated by the AWD and CMI, respectively [46]. It was highlighted that the CMI captured fewer drought events than the AWD. This fact can be linked to the unrealistic zero or near zero values from the CMI under dry periods, as previously discussed.

The SWDI_Rawls captured between 88.9% and 100% of the drought weeks, and the SWDI_CYL captured between 78.6% and 100%, except in SA101, SA102, SG02 and Salamanca (Figure 7a,b). These results were higher than the 67.3–78.9% obtained by the SWDI computed from the root-zone SM stations of REMEDHUS (2009–2013) and the estimated FC and WP from the water retention curves [33]. Similarly, they were also higher than the 59.6–71.2% obtained by the SWDI computed from the SMOS L2 SM data and the specific PTFs derived over REMEDHUS (June 2010–December 2014) [46]. This may be due to the larger region analyzed in this study, which includes a higher variability of soil textures than the previous research about SWDI in REMEDHUS. Both SWDI approaches captured more drought weeks than the AWD and the CMI, except in SA101, SA102, SG02 and Salamanca. At these stations, the number of drought weeks identified by both SWDI approaches (from 63.6% to 81.8% for SWDI_Rawls and from 27.9% to 68.0% for SWDI_CYL) was lower than at the other stations. This was explained by the very coarse soil texture of SA101, SA102, SG02 and Salamanca, which resulted in a low estimation of the FC value and, consequently, a high SWDI compared with the other stations (see Figure 1 and Figure 2).

The number of drought weeks captured by both SWDI approaches was always higher than the number captured by the AWD, CMI, SMADI, SMDI or SWetDI. This could be related to the dry bias of the SMOS SM products, which was found in several validation exercises using a previous version of the SMOS L2 SM data (v.551), and the subsequent SMOS BEC L3 SM (v.1) and L4 SM (v.2 and v.3) products in several studies over REMEDHUS and Castilla y León regions [59,64]. Although this study used a new version of the data (v.620), this dry bias is still present in some regions [37]. The SMOS dry bias has a direct influence on the SWDI. In contrast, the SMADI uses normalized SM data, and the SMDI and SWetDI use an SM anomaly with the bias removed. As previously stated, the SMOS is only capable of retrieving the SM content of the soil top layer, which surely underestimates the root-zone SM. In the near future, the availability of root-zone SM estimates [92,93] is expected to improve the current estimation based on surface SM data for agricultural drought monitoring. Nevertheless, recent research has shown the capability of using satellite surface SM estimations for drought monitoring, e.g., using the surface SM observations from the Advanced Microwave Scanning Radiometer for Earth (AMSR-E) observation system and SMAP, after applying a cumulative distribution function (CDF) matching in situ surface SM data to remove the systematic bias and dynamic range differences [79,94].

The SMADI captured from 27.7% to 69.3% of the drought weeks (Figure 7a,b). These results were lower than those obtained by the AWD and both SWDI approaches but were similar to those obtained by CMI, SMDI and SWetDI. A similar number of drought weeks was identified by the SMDI and SWetDI, from 44.2% to 59.1% and from 41.7% to 56.5%, respectively (Figure 7a,b). This could be explained because both indices use SM or SWetI anomalies, as previously discussed. These results were also lower than those obtained by AWD and both SWDI approaches.

3.4. Categorical Statistical Analysis

The categorical statistics of the drought indices are shown in Figure 8. First, it is important to mention that even though the AWD and CMI were considered the reference dataset, they should not be considered as the absolute “ground-truth”. They were assumed to be benchmarks because they are a good proxy of the soil water deficit (the AWD) and a widely used agricultural drought index (the CMI). Therefore, the results inferred from the statistical analysis only showed the level of similarity of the drought indices with respect to the AWD and CMI.

Considering the POD (Figure 8a,b), both SWDI approaches estimated a higher fraction of the AWD and CMI drought events than the other indices (median values of 0.99 and 1 for SWDI_Rawls, 0.96 and 0.98 for SWDI_CYL, 0.57 and 0.78 for SMADI, 0.54 and 0.67 for SMDI, and 0.55 and 0.58 for SWetDI for AWD and CMI, respectively). The SWDI results were in line with those of a recent study [78], where a POD value of 0.98 was observed between SWDI and AWD. When analyzing the POFD (Figure 8c,d), the SMADI displayed the lowest values (median values of 0.10 and 0.19 for the AWD and the CMI, respectively), followed by the SMDI (0.35 and 0.36), the SWetDI (0.43 in both cases), the SWDI_CYL (0.76 and 0.88) and the SWDI_Rawls (0.92 and 0.97).

The FAR of the AWD was lower than 0.20 for all indices (medians of 0.18 for SWDI_Rawls, 0.16 for SWDI_CYL, 0.04 for SMADI, 0.13 for SMDI and 0.15 for SWetDI, Figure 8e), but the FAR of the CMI was higher (medians of 0.55 for SWDI_Rawls, 0.54 for SWDI_CYL, 0.23 for SMADI, 0.41 for SMDI and 0.50 for SWetDI, Figure 8f). When analyzing the FB, both SWDI approaches slightly overestimated, while the other indices underestimated the number of drought weeks with respect to the AWD (medians of 1.21 for SWDI_Rawls, 1.16 for SWDI_CYL, 0.59 for SMADI, 0.61 for SMDI and 0.64 for SWetDI, Figure 8g). By contrast, all indices overestimated the number of drought weeks with respect to the CMI (medians of 2.25 for SWDI_Rawls, 2.15 for SWDI_CYL, 1.03 for SMADI, 1.11 for SMDI and 1.15 for SWetDI, Figure 8h). The FAR and the FB results agreed with those previously shown in Figure 7.

Both the POD and POFD referred to the AWD and the CMI were very similar, whereas the FAR and the FB exhibited higher values when the CMI was assumed as the reference, which is in line with results shown in Figure 7. Considering all categorical statistics, both SWDI approaches were similar to the AWD. In contrast, the SMADI, the SMDI and the SWetDI were similar to the CMI.

3.5. Spatial Comparison

The SMOS BEC L4 SM and MODIS SWetI maps over the Castilla y León region at the two selected weeks are shown in Figure 9. The week of 8–14 April 2012 was drier than 9–15 April 2013 (Figure 8a,b), with values from 0.10 to 0.20 m³/m³. In addition, some parts of the northwest and the center were drier than others (Figure 9a). These spatial features did not appear in the SWetI map of the same week (Figure 9c). Although the purpose of the SWetI is to mimic the soil water content behavior, the SWetI did not seem to reflect the same pattern as the SM in this week. The SM map of 9–15 April 2013 exhibited a high SM content and high variability (from 0.12 to 0.40 m³/m³) and showed that the southeast part was wetter than the rest of the area (Figure 9b). It was clearly observed that the spatial patterns of the SM during 9–15 April 2013 were similar to those of the SWetI of the same week (Figure 9d). In this case, the SWetI and the SM maps seemed to capture similar information about the soil water content.

The maps of the agricultural drought indices over the Castilla y León region during the two selected weeks are shown in Figure 10. To show a homogeneous legend, the drought intensity is depicted with different colors according to the thresholds of each index. As mentioned earlier, a very intense meteorological drought occurred in 2012 on the Iberian Peninsula [83].

During the week of 8–14 April 2012, the SWDI_Rawls displayed severe and extreme drought in the north, whereas the south exhibited no drought to moderate drought conditions, except for small areas of the west and south that showed severe and extreme drought conditions (Figure 10a). In the case of SWDI_CYL, the spatial patterns were similar to those found with the SWDI_Rawls in the north. However, no drought was estimated by the SWDI_CYL in most areas of the south (Figure 10c). These results agreed with the spatial distribution of SM in the same week, as observed in Figure 9a. The difference between the north and the south was related to the different soil textures of these areas, with lower sand content in the north than in the south, as observed in Figure 2. Therefore, the SWDI seemed to be sensitive to the different soil types, which agrees with its rationale and is in accordance with the definition of agricultural drought. The low values of both SWDI approaches in the south may be due to the different water holding capacities of these soils, which can be inferred from the soil characteristics of this part of the study area (Figure 2).

The SMADI exhibited few scattered pixels with mild drought conditions in most areas during 8–14 April 2012 (Figure 10e). The drought intensity indicated by the SMADI was lower than that displayed by both SWDI approaches. This could be because the MTCI/VCI term of the SMADI modulates the soil water content through temperature and vegetation conditions. During this week, there were not extreme LST values, as occurs during the summer seasons in this region. This fact may reduce the SMADI drought intensity. In addition, there were some small zones where the SMADI estimated severe and extreme drought, which corresponded to areas of irrigated crops. At this stage of the growing season, these small zones only had bare soil. The absence of vegetation is characterized by a very low NDVI, which could, in turn, produce very high SMADI values. In this case, these high SMADI values were unrealistic and not associated with drought. This was because some factors (e.g., fire, flooding, hail, pests, plant disease or land cover change) can produce a vegetation signal that mimics drought conditions [95].

During 8–14 April 2012, the SMDI showed extreme drought in the northwest, severe drought in some parts of the north, and a varying pattern ranging from no drought to moderate drought in other areas (Figure 10g). A similarity was observed between the spatial patterns of the SMDI and both SWDI approaches. Nevertheless, the SMDI presented a lower drought intensity than the SWDI and a higher drought intensity than the SMADI. The SMDI computation employed the weekly SM anomaly data instead of the absolute SM information, and the variability of the SMDI was very high from one week to the next. Due to that, the similarity between the SMDI and the SM maps during this particular week (Figure 8a and Figure 9g) could not be seen in other weeks, as occurs during 9–15 April 2013 (Figure 8b and Figure 9h).

The SWetDI showed a range from no drought conditions to mild drought in all areas, and only a few pixels with moderate drought occurred for the week of 8–14 April 2012 (Figure 10i). The SWetDI exhibited a lower drought intensity than the SMDI and both SWDI approaches and a higher drought intensity than the SMADI. Although the SWetDI used only a monthly anomaly, its spatial patterns were very similar to those obtained with the SWetI, as shown in Figure 9b. However, the SWetI estimation did not always capture the same soil water content information as the SM estimation, as previously discussed.

During the week of 9–13 April 2013, the difference between the north and south areas was also clearly observed with both SWDI approaches (Figure 10b,d). In addition, there were some areas with severe and extreme drought in the north in both SWDI approaches. This could be related to the fine texture of the soil in the north part of the study region. Although fine soils have more soil surface water content, this water would not have been available to plants. In general, the SWDI_Rawls estimates a higher drought intensity than the SWDI_CYL. The SMADI considered no drought in all areas of the study region (Figure 9f). In this case, the unrealistically high SMADI values in bare soil zones did not appear, probably because of the higher SM contents of this week compared that of 8–14 April 2012. No drought was also estimated by the SMDI over the entire region during 9–13 April 2013 (Figure 10f). In addition, the spatial patterns of the SMDI were not similar to those of the SM for this particular week, as observed in Figure 8b. The SWetDI displayed no drought in the majority of the region, but there were some areas of mild drought, especially in the northeast (Figure 10j).

4. Conclusions

During the last decade, the great importance of SM in agricultural drought monitoring led to the development of a variety of agricultural drought indices that mainly use this variable as drought indicator under different points of view. In this study, the performance of four innovative satellite SM-based agricultural drought indices was evaluated over a large agricultural area in northwest Spain. The indices were weekly estimated from June 2010 to December 2016 using the SMOS BEC L4 SM and/or the MODIS SR and LST data at 1 km spatial resolution. In the case of the SWDI, two approaches (SWDI_Rawls and SWDI_CYL) that were obtained from different PTFs were investigated separately. A temporal analysis was conducted by comparing these indices with two well-known agricultural drought indices: the AWD and the CMI. In addition, the spatial patterns of the estimated drought index maps were also compared under dry and wet SM conditions at the beginning of the growing season.

The SWDI_Rawls, SWDI_CYL and SMADI adequately reproduced the seasonality of the AWD and CMI. Conversely, this seasonality did not appear in the SMDI and SWetDI time series. On the one hand, the length of the available SM data, or its surrogate, was considered in all SM-based drought indices computed, except in the SWDI. This issue had a noticeable influence on the SMDI computation. On the other hand, the size of the study area influenced the SWetI estimation and, consequently, the SWetDI computation. Nevertheless, the SMDI trend was capable of capturing the very intense meteorological drought of 2012, which was not captured by the SWetDI trend.

In general, the SWDI_Rawls worked equally well as the SWDI_CYL. Therefore, the SWDI could be calculated in any region with available soil databases by applying suitable standards and adequately tested PTFs. However, both SWDI approaches were sensitive to the accuracy of the FC and WP estimations, which, in turn, were sensitive to the robustness of the soil database, especially when there was a very low SM content. By contrast, the SMADI was not directly affected by very low SM, because it includes additional LST and NDVI information.

In comparison with the AWD, both SWDI approaches obtained the highest correlation (R ≈ 0.68–0.80), followed by SMADI (R ≈ −0.60 to 0.83). At the same time, very low and non-significant correlations were found with the SMDI (R ≤ 0.3) and the SWetDI, respectively. Similar results were obtained when the drought indices were compared with the CMI (R ≈ 0.48–0.70 for both SWDI approaches, R ≈ −0.32 to −0.69 for SMADI, R ≈ 0.39 to 0.44 for SMDI and R ≤ 0.28 for SWetDI). Both SWDI approaches overestimated the drought weeks compared to the AWD and the CMI. The SMADI, the SMDI and the SWetDI underestimated the drought weeks compared to the AWD, whereas they overestimated the drought weeks compared to the CMI. In terms of similarity, the SWDI agreed with the AWD, and the SMADI, SMDI and SWetDI agreed with the CMI.

When analyzing the estimated drought index maps, a spatial pattern with two distinctive soil property areas was found in the SWDI maps, which displayed severe and extreme drought in fine-textured soils. The SMADI only estimated a few scattered pixels with mild drought under very low SM conditions, since its intensity was modulated by the temperature and vegetation conditions. No distinctive patterns were observed in the SMDI and the SWetDI maps.

The estimation of the soil water content from the optical VIS/IR observations through the LST-NDVI triangle does not always capture the same information as the passive microwave SM estimations. Then, the use of the SM variable itself is recommended. In addition, the absolute SM is preferred as agricultural drought indicator over the SM anomalies, since SM anomalies removed the inherent seasonality of the soil water content. Therefore, the SWDI and the SMADI would be good options for developing a feasible tool dedicated to water resource management and drought monitoring in many cropland regions worldwide.