Spatiotemporal Analysis of Rainfall and Droughts in a Semiarid Basin of Brazil: Land Use and Land Cover Dynamics

Abstract

Precipitation estimation is a challenging task, especially in regions where its spatial distribution is irregular and highly variable. This study evaluated the spatial distribution of annual rainfall in a semiarid Brazilian basin under different regimes and its impact on land use and land cover dynamics. Climate Hazards Group InfraRed Precipitation with Stations (CHIRPS) records and observed data from 40 weather stations over a time series of 55 years were used, in addition to the Standardized Precipitation Index. Spatiotemporal analysis was carried out based on geostatistics. Remote sensing images were also interpreted for different rainfall regimes using the Normalized Difference Vegetation Index and Modified Normalized Difference Water Index. The Gaussian semivariogram model best represented the rainfall spatial structure, showing strong spatial dependence. Results indicated that rainfall amount in the basin significantly increases with elevation, as expected. There is high variation in the dynamics of water storage that can threaten water security in the region. Our findings point out that the application of geostatistics for mapping both the annual precipitation and the Standardized Precipitation Index provides a powerful framework to support hydrological analysis, as well as land use and land cover management in semiarid regions.

1. Introduction

Precipitation is an essential factor in hydrological processes (e.g., erosion, flooding, and aquifer recharge) [1,2,3,4]. In that regard, the spatiotemporal distribution of rainfall has a major impact on the hydrological cycle, water security, agriculture, water resource management, and on several economic and productive activities [5,6,7].

Despite the relevance of the rainfall distribution in water availability, reliable and high-quality long-term rainfall time series remain scarce in South America, particularly in northeastern Brazil (NEB) [8]. Moreover, precipitation estimation is particularly challenging in the Brazilian semiarid region due to its highly variable spatial distribution and irregular periodicity, which has caused frequent droughts over the years [9,10]. Thus, applying traditional and widely adopted meteorological indices for monitoring droughts in such areas, such as the Standardized Precipitation Index (SPI), is usually a hard task that requires data from several rainfall stations [11,12,13,14]. The usefulness of SPI by analyzing long-term rainfall data from several drainage basins in the eastern Mediterranean region was evaluated by Tsesmelis et al. [15], emphasizing the importance of SPI as a tool for assessing changes in drought and wetness episodes, enabling the identification of areas requiring water management strategies to improve water security and to reduce the impact of hydrological extremes.

Remote sensing is an important tool for hydrological research, allowing for the estimation of rainfall using satellite images with high spatiotemporal resolution. The Tropical Rainfall Measuring Mission (TRMM), Climate Hazards Group InfraRed Precipitation with Stations (CHIRPS), and Precipitation Estimation from Remotely Sensed Information using Artificial Neural Networks (PERSIANN) have been widely adopted for rainfall analysis on a basin scale [8,10,16,17]. The CHIRPS product offers a high spatial resolution of approximately 5.3 km and a long temporal period (from 1981) for validation using in situ rain gauge measurements [8,16,18]. In addition, remote sensing can also be used to monitor landscape and biophysical changes in semiarid regions through the evaluation of land use and land cover (LULC) map spatial association with vegetation indices, such as the Normalized Difference Vegetation Index (NDVI) and the Modified Normalized Difference Water Index (MNDWI). These indexes allow the estimation of land use and cover, enabling the characterization of spatiotemporal dynamics and quantification of vegetation conditions and water features [5,19,20,21,22,23,24]. Additionally, the NDVI has proved to be a valuable tool in the analysis of irrigation in agricultural areas, allowing farmers and researchers to assess soil salinity and monitor crop growth and water stress to optimize yield and conserve resources [25,26]. Due to the high resolution and wide availability of remote sensing products for hydrological studies, there has been a surge of interest from the scientific community to validate such indirect information for different basins around the world [18,24,27]. Moreover, remote sensing products are highly valuable for climate change investigations, and especially for devising adaptive strategies to promote sustainability, particularly in semiarid regions.

The semiarid region of Brazil is known for its highly unpredictable weather patterns, featuring short wet periods and long dry spells. This region is considered the wettest semiarid region worldwide, as evidenced by Marengo et al. [28]. Furthermore, Montenegro and Ragab [29] demonstrated that rivers in this region predominantly have non-perennial hydrological regimes. One of the most prominent rivers in the region is the São Francisco River, with its tributary, the Brígida River, being the focus of this study. The Caatinga vegetation biome dominates the region, which is characterized by deciduous cover that sheds part of its canopy during the dry season.

Typically, rainfall data estimates are obtained by using regression techniques based on time series analysis [30]. The application of geostatistical modeling through kriging takes into account the structure of spatial variation and usually provides high-precision maps, which are important for the study of extreme events [27,31,32,33]. This method is especially important for the NEB, where the density of rainfall networks is usually low, and proper characterization of rainfall spatial variability is limited due to the lack of observational data [34]. Among the geostatistical estimation methods, kriging is often used in the literature [32,35,36,37].

The spatiotemporal variability of monthly rainfall data at 269 rainfall stations from 1994 to 2014 in the Paraíba State, located in the NEB, was analyzed by Medeiros et al. [37], who found that 85% of the total rainfall variability could be well represented using kriging, providing a reliable spatial distribution of regions where rainfall conditions varied from scarcity to excess. Regarding the spatiotemporal variability of annual rainfall, Silva et al. [35] found substantial annual spatiotemporal variability and a negative trend (decrease) for most rainy seasons in the Epitácio Pessoa reservoir watershed, Paraíba State, during the study period (1963–1991) using the kriging technique. Additionally, for the Brazilian semiarid region, Araújo et al. [27] investigated the spatiotemporal distribution of precipitation in the Brígida River basin, Pernambuco State, using kriging, providing an adequate precipitation estimate and classification of the aridity condition; however, the availability of observational data was a limiting factor. Furthermore, due to the high spatial variability of precipitation patterns related to topographic and climatic factors [38,39], a study conducted by Adhikary et al. [39] compared three geostatistical methods, finding that the cokriging method with elevation as a secondary variable was significant in Australian watersheds with complex mountainous terrain. Similarly, relationships between elevation and rainfall have been poorly studied in the semiarid regions of Brazil.

Despite the effectiveness of kriging in reducing the variance of the estimation error, it is strongly dependent on the quantity and quality of data, as well as their spatial distribution. However, it may not always accurately reproduce the existing spatial correlation, which can result in unrealistic spatial distributions [40,41]. An alternative estimation method is based on a stochastic simulation approach, which, by definition, reproduces the statistical behavior of the phenomenon. The most widely used simulation technique for forecasting and uncertainty assessment is the sequential Gaussian simulation (SGS) method [42]. The main advantage of the SGS is to produce alternative and equally probable stochastic realizations that reproduce model statistics rather than reduce local forecast variation [43]. In a study of the spatial variability of extreme rainfall in Mexico, Méndez-Venegas et al. [41] used stochastic simulation techniques with a geostatistical approach, highlighting the Gaussian simulation as a powerful tool for studying rainfall variability and modeling its spatial distribution. Another method that seeks to improve the accuracy of spatial estimates is cokriging, which typically reduces the forecast error variance and outperforms kriging when the chosen secondary variable is highly correlated with the primary variable [44,45].

Although a considerable number of studies have addressed the spatial variability of rainfall patterns in the literature, there is still a lack of research focusing on the spatial dependence structures for various hydrological regimes (from dry to wet years) based on the Standardized Precipitation Index (SPI) and drought indexes, especially in the Brazilian semiarid region. Furthermore, there is a need for a joint analysis of the relationship between water and vegetation in the region. Indeed, the SPI values for a 12-month period are strongly linked to streamflow regimes, reservoir levels, water security, and long-term management [12,13,34,46].

Therefore, the aim of this research was to analyze the distribution of rainfall in a basin area of the Brazilian semiarid region, using both observational and remote sensing data. The study seeks to provide a methodological framework for characterizing the spatial scales of rainfall dependence in long-term semiarid regions. This study applies both deterministic and stochastic geostatistical procedures to map and examine the correlation of rainfall with elevation, vegetation cover, and water availability in order to assist decision-making related to water security, vegetation–water nexus studies, reforestation, climate change initiatives, and public policies. Additionally, it assesses the impacts of rainfall variability on land use and land cover (LULC) dynamics using remote sensing data and vegetation indices.

2. Materials and Methods

2.1. Characterization of the Study Site

The study was conducted in the Brígida River basin (Figure 1), located in the High Sertão mesoregion of the Pernambuco State, Brazil (07°19′02″ to 08°36′32″S and 39°17′33″ to 40°43′06″W), with an average elevation of 557 m. The basin area belongs to the semiarid region of Pernambuco, covering 13,495.73 km², which corresponds to 13.73% of the total surface of the state. It comprises 15 municipalities and is part of the São Francisco River basin area [47]. The study basin is one of the tributaries of the São Francisco River Integration Project (PISF), which plays a strategic role in ensuring future water availability in the semiarid northeast region. The PISF main objective is to provide water supply for the socioeconomic development of the most drought-vulnerable northeast states of Brazil (Ceará, Paraíba, Pernambuco, and Rio Grande do Norte) [48].

According to the Köppen–Geiger climate classification, the region is classified as type BSh (hot semiarid climates) and Aw (tropical with dry winter) [49,50], as shown in Figure 1. The monthly mean temperature is 26.1 °C, indicating an average potential evapotranspiration of 175.1 mm per month and an annual cumulative value greater than 2000 mm. The mean total annual precipitation is 608.0 mm (29 years; from 1991 to 2020), with the wet season occurring from December to April (with an average total rainfall of 558.5 mm) [51]. The predominant vegetation in the basin is the hyperxerophilous Caatinga, and the predominant soils are Oxisols, Ultisols, Alfisols, Entisols, and Aridisols [52].

2.2. Rainfall Data

Data from the Climate Hazards Group InfraRed Precipitation with Station (CHIRPS) was used from 1981 to 2017 (https://developers.google.com/earth-engine/datasets/catalog/UCSB-CHG_CHIRPS_DAILY, accessed on 9 May 2023), which includes data from 50°S to 50°N latitudes. CHIRPS combines in situ climatology and satellite imagery with a spatial resolution of 0.05° (about 5.3 km) spanning 35 years [8,9,10,11,12,13,14,16,53], and has already been validated for the Brazilian semiarid region [8]. The CHIRPS dataset was validated for the selected basin using observed data from 40 rainfall stations, with an annual time series of 55 years (1962 to 2017), provided by the Pernambuco State Water and Climate Agency (APAC). These data were tested for homogeneity and missing data by APAC [54].

To evaluate the performance of CHIRPS against the observed precipitation data, a pixel-by-pixel comparison was conducted based on the location of the rainfall stations, and various statistical indices were used, including the coefficient of determination (R²), Pearson’s correlation (r), Nash–Sutcliffe Efficiency Coefficient (NSE), and percentage of bias (PBIAS). CHIRPS data, with its high resolution, has the potential to provide a better representation of local-scale spatial variation of rainfall in the basin.

2.3. Standard Precipitation Index (SPI)

The Standard Precipitation Index (SPI) [11] was employed to evaluate the temporal variability of cumulative annual precipitation (SPI-12), commonly adopted to monitor drought and excessive rainfall conditions [12,13,14]. We used the time series of precipitation observed at rainfall stations from 1962 to 2017, which was later supplemented with CHIRPS data from 1981 to 2017. The SPI is given in thresholds, which correspond to the number of standard deviations that the observed cumulative precipitation deviates from the climatological mean. The precipitation classification used was described by McKee et al. [11] and grouped according to Table 1, which ranges from the most extreme deficit to excess precipitation situations. We selected years for each precipitation class from 1981 onwards, based on their classification using the SPI-12, to carry out the statistical analysis.

2.4. Statistical Analysis

The CHIRPS data used, covering the entire basin area, consisted of 930 pixels, and were collected from 1981 to 2017. These data were subjected to classical statistics to verify their behavior in terms of statistical measures of position and variability. According to the coefficient of variation (CV) values, variability was classified as low (CV ≤ 12%), intermediate (12% < CV ≤ 60%), and high variability (CV > 60%), following Warrick and Nielsen [55]. Additionally, the data distribution was evaluated in relation to the normal distribution using the Kolmogorov–Smirnov test at a 5% significance level. The data were then used to estimate the average semivariogram of the study area. Trend analysis was performed by estimating the trend surface using a quadratic polynomial function and examining the coefficients of determination obtained. Trends were considered to exist for R² values greater than or equal to 0.7. In such cases, residual values were adopted for the analysis.

Geostatistical analysis was conducted using GS+ software. The sampling locations corresponded to the center point of the CHIRPS pixels and were supplemented with observed values from stations. Spatial dependence was evaluated using the classic semivariogram built from the estimate of classical semivariance, as given by Equation (1) [56].

(1)$\mathsf{\gamma }\left(\mathrm{h}\right)=\frac{1}{2\mathrm{N}\left(\mathrm{h}\right)}{\displaystyle \sum }_{\mathrm{i}=1}^{\mathrm{N}\left(\mathrm{h}\right)}{[\mathrm{Z}\left({\mathrm{X}}_{\mathrm{i}}\right)-\mathrm{Z}\left({\mathrm{X}}_{\mathrm{i}}+\mathrm{h}\right)]}^2$

where γ(h) is the estimated value of the semivariance of the experimental data; Z(X_i + h) and Z(X_i) are the observed values of the regionalized variable; N(h) denotes the number of pairs of measured values separated by distance h; and X_i and X_i + h refer to two data points (locations where the variable Z is measured) that are a lag distance, h, apart [57].

According to Landim [57], the sums required to calculate γ(h) must include a sufficient number of pairs to produce consistent results, with a minimum number of 30 pairs per distance class adopted [58]. After evaluating the experimental semivariogram, exponential, spherical, and Gaussian theoretical models were tested, and the parameters C₀ (nugget effect), C₀ + C₁ (sill), and A (range of spatial dependence) were estimated. The three theoretical models considered are presented below in Equations (2)–(4) [40]:

• Exponential Model:

(2)$\mathsf{\gamma }\left(\mathrm{h}\right){=\mathrm{C}}_0{+\mathrm{C}}_1\left\{1-\exp \left[-\left(\frac{\mathrm{h}}{\mathrm{A}}\right)\right]\right\}, \mathrm{h}\ne 0$

• Spherical Model:

(3)$$\begin{gathered} \mathsf{\gamma }\left(\mathrm{h}\right){=\mathrm{C}}_0{+\mathrm{C}}_1, \mathrm{h}\ge \mathrm{A} \\ \mathsf{\gamma }\left(\mathrm{h}\right){=\mathrm{C}}_0{+\mathrm{C}}_1\left[1.5\left(\frac{\mathrm{h}}{\mathrm{A}}\right)-0.5{\left(\frac{\mathrm{h}}{\mathrm{A}}\right)}^3\right], 0<\mathrm{h}<\mathrm{A} \end{gathered}$$

• Gaussian Model:

(4)$\mathsf{\gamma }\left(\mathrm{h}\right){=\mathrm{C}}_0{+\mathrm{C}}_1\left\{1-\exp \left[-{\left(\frac{\mathrm{h}}{\mathrm{A}}\right)}^2\right]\right\}, \mathrm{h}\ne 0$

where C₀ is the nugget effect; C = C₀ + C₁ is the sill; A is the range; and h denotes the distance between points.

The degree of spatial dependence (DSD) was assessed based on the methodology proposed by Cambardella et al. [59], which uses the relationship between the nugget effect and the sill of the adjusted semivariogram, classifying the spatial dependence as strong, moderate, or weak, as observed in Equation (5). Values below 25% are characterized as strong spatial dependence, values between 25% and 75% are considered moderate, while values above 75% indicate weak spatial dependence.

(5)$\mathrm{DSD}(\%)=\frac{{\mathrm{C}}_0}{{\mathrm{C}}_0{+\mathrm{C}}_1}$

The best models were selected by analyzing the coefficient of determination (R²) and performing cross-validation using the Jack-Knifing criterion [60], assuming a mean error close to zero and a standard deviation close to one.

2.4.1. Universal Kriging (UK)

After modeling the semivariograms, the values were interpolated at non-sampled sites using Universal Kriging (UK). This technique uses both correlations with auxiliary maps and spatial correlation, where the drift (or trend) is modeled as a function of coordinates [61], and then kriging is performed on the local average residuals [39].

2.4.2. Sequential Gaussian Simulation (SGS)

The sequential Gaussian simulation (SGS) was applied to different classifications of SPI. At each location, the SGS algorithm randomly selects a simulated value from an estimated conditional cumulative distribution function. The distribution function is determined by the mean and variance of the kriging calculated from neighborhood information. Running the SGS requires transforming the original data into a Gaussian distribution, which is accomplished by transforming the normal score. The cumulative density function is characterized by the mean and covariance [62], assuming a random Gaussian field. By repeating these sequential steps with different random paths, the SGS algorithm can provide various realizations of the spatial distribution of the simulated values. In this study, to obtain an accurate probability calculation, the SGS algorithm was run with ten, one hundred, and one thousand runs using GS+ software.

2.4.3. Ordinary Cokriging

After modeling the semivariograms, the values were interpolated at unsampled locations to estimate the primary variable using the ordinary cokriging (CKO) procedure, which is a modification of the simple kriging method. CKO offers the advantage of using more than one variable in the estimation process, thereby improving the estimation of the primary variable by assuming correlation between variables [44]. In this study, rainfall and elevation were considered the primary and secondary variables, respectively, in the CKO method. Thus, CKO aimed to establish a suitable model for cross-continuity and spatial dependence between the primary (rainfall) and secondary (elevation) variables. This correlation between variables is known as interregionalization [43], which can be quantified by cross semivariogram or cross-covariance, as seen in Equation (6).

(6)$\mathsf{\gamma }\left(\mathrm{h}\right)=\frac{1}{2\mathrm{N}\left(\mathrm{h}\right)}{\displaystyle \sum }_{\mathrm{i}=1}^{\mathrm{N}\left(\mathrm{h}\right)}{[\mathrm{Z}}_1\left({\mathrm{X}}_{\mathrm{i}}\right)-{\mathrm{Z}}_1\left({\mathrm{X}}_{\mathrm{i}}+\mathrm{h}\right){\left]\right[\mathrm{Z}}_2\left({\mathrm{X}}_{\mathrm{i}}\right)-{\mathrm{Z}}_2\left({\mathrm{X}}_{\mathrm{i}}+\mathrm{h}\right)]$

where γ(h) is the estimated value of the cross semivariogram of the experimental data; Z(X_i + h) and Z(X_i) are the observed values of the regionalized variable; Z₁ and Z₂ are spatially correlated variables; N(h) is the number of pairs of measured values separated by a distance, h; and X_i and X_i + h refer to two data points (locations where the variable Z is measured) that are a lag distance, h, apart [57].

The elevation was obtained through a Digital Elevation Model (DEM) using the Shuttle Radar Topography Mission (SRTM) product, which was acquired from the Google Earth Engine (GEE) with a spatial resolution of 30 m.

2.5. Assessment of Geostatistics Methods

In this study, the performance of all geostatistical methods for rainfall estimation is compared based on the values of mean bias error (MBE), root mean squared error (RMSE), and coefficient of determination (R²) between observed and estimated rainfall values, according to Adhikary et al. and Hengl et al. [39,61].

The geostatistical method with the lowest MBE and RMSE values and the highest R² value is selected as the best geostatistical method.

2.6. Analysis of Biophysical Indexes and Land Use and Land Cover (LULC)

To characterize vegetation cover in the basin and highlight water features in a semiarid region for different years, the Normalized Difference Vegetation Index (NDVI) and Modified Normalized Difference Water Index (MNDWI) were calculated. Images from Landsat 5 TM and Landsat 8 OLI sensors were used, capturing different rainfall regimes (MD, N, MW, and EW) during selected years after the month of September, as detailed in Table 2. Images with less than 10% cloud coverage were processed on the Google Earth Engine (GEE) platform (https://developers.google.com/earth-engine/datasets/catalog/landsat, accessed on 9 May 2023), provided by the National Aeronautics and Space Administration (NASA). The images were then transferred to Google Earth Engine (GEE) where they were atmospherically corrected for surface reflectance. All images were processed using GEE, classified, and mapped using QGIS software version 3.22.

The Normalized Difference Vegetation Index (NDVI) is commonly used and disseminated in the literature. It is calculated as the ratio between the difference of the near-infrared and red bands’ reflectance (ρ_NIR and ρ_R, respectively) and their sum, according to Equation (7) [19].

(7)$\mathrm{NDVI}=\frac{{\mathsf{\rho }}_{\mathrm{NIR}}-{\mathsf{\rho }}_R}{{\mathsf{\rho }}_{\mathrm{NIR}}+{\mathsf{\rho }}_R}$

In turn, the Modified Normalized Water Difference Index (MNDWI) proposed by Xu [64] aims to capture water bodies on the earth’s surface and identify water courses. The MNDWI can be obtained by calculating the ratio between the difference in reflectance of the green (ρ_G) and mid-infrared (ρ_MIR) bands, divided by their sum (Equation (8)).

(8)$\mathrm{MNDWI}=\frac{{\mathsf{\rho }}_{\mathrm{G}}-{\mathsf{\rho }}_{\mathrm{MIR}}}{{\mathsf{\rho }}_{\mathrm{G}}{+\mathsf{\rho }}_{\mathrm{MIR}}}$

This index minimizes information about vegetation and accurately highlights water bodies such as rivers and reservoirs, where positive MNDWI values indicate pixels filled with water [64]. An estimate was made of the volume of water stored in the surface reservoirs of the basin using Equation (9), proposed by Molle [65], who studied more than 400 reservoirs in the Brazilian semiarid region and proposed a model that considers only the area of formed water bodies. The proposed model is highly accurate, with a coefficient of determination of 0.946.

(9)${\mathrm{V}}_{\mathrm{R}}=1{.98 \mathrm{S}}_{\mathrm{x}}{}^{1.001}$

where V_R is the volume of the reservoir (m³), and S_x is the area of the water bodies of the reservoir (m²), as adopted by Silva et al. [24].

The land use and land cover (LULC) maps were obtained using the MapBiomas Brazil platform (https://mapbiomas.org/, accessed on 9 May 2023) for the same years analyzed in the biophysical indices (Table 2) to characterize the study area and investigate the impacts of land use changes and climate variability. The extracted data were from the latest available collection, collection 7.1, with a spatial resolution of 30 m and a catalog with 27 classes, in which the following major classes were prioritized: Forest, Non-forest natural formation, Farming, Non-vegetated area, and Water [66]. Additionally, the area values (km²) corresponding to each of the major classes studied were extracted to quantify the evolution of these areas over the different rainfall regimes.

2.7. Flowchart of the Data Analysis and Processing Steps

In order to represent the stages of statistical analysis, geostatistical modeling, development of kriging, sequential Gaussian simulation and cokriging maps, and analysis of biophysical indices, Figure 2 shows the flowchart of the processing steps.

3. Results and Discussion

3.1. Rainfall via CHIRPS and SPI

Figure 3 depicts the correlation and accumulated density curve of observed rainfall data and CHIRPS precipitation estimates. The percentage of trend (PBIAS) was 9.2%, indicating a “very good” performance (<10%). The coefficient of determination exhibited a high-quality fit (R² = 0.938). The Nash–Sutcliffe Efficiency Coefficient (NSE) was 0.63, classified as “satisfactory” (˃0.50), and the Pearson correlation coefficient (r) was 0.82, further demonstrating the correlation of CHIRPS with weather stations. Additionally, the angular parameter of the linear regression model between observed and estimated rainfall data is very close to one, meaning no severe underestimation or overestimation when utilizing the estimated data to represent rainfall dynamics in the basin. These results align with those obtained by Paredes-Trejo et al. [8], where they validated satellite rainfall estimates based on CHIRPS for northeast Brazil via automatic meteorological stations and achieved a Pearson correlation coefficient (r) of 0.941, PBIAS = −3.58%, and NSE = 0.886.

Based on the observed rainfall data from stations before 1981 and after 1981, along with CHIRPS data on total annual rainfall, the Standard Precipitation Index (SPI-12) was calculated. The SPI-12 classifies years into Extremely Dry (ED), Severely Dry (SD), Moderately Dry (MD), Normal (N), Moderately Wet (MW), Severely Wet (SW), and Extremely Wet (EW) categories (Table 3). The majority of the years classified as N years; followed by MD years; MW years; ED years; SW years; and, with only one occurrence each, SD years (1983) and EW years (1985).

3.2. Statistical Analysis of Annual Rainfall

Table 4 presents the descriptive statistics of the annual rainfall data estimated by CHIRPS for the years 1993 (ED), 1983 (SD), 2017 (MD), 2014 (N), 2004 (MW), and 1985 (EW) in 930 pixels. The coefficient of variation (CV) for these years was intermediate (12% < CV ≤ 60%), based on the SPI classification in Table 4, according to Warrick and Nielsen [55]. Additionally, Medauar et al. [36] evaluated the spatial behavior and temporal stability of the mean monthly rainfall in the Bahia State, Brazil, which is part of the NEB, from 1975 to 2011, and found intermediate CV for all analyzed months, which supports the findings of this study.

Figure 4 presents boxplots for annual rainfall estimated by CHIRPS for the years 1993 (ED), 1983 (SD), 2017 (MD), 2014 (N), 2004 (MW), and 1985 (EW). Rainy years (1985 and 2004) showed greater data dispersion, unlike dry years such as 1983, which presented lower dispersion, corroborating the findings of Medeiros et al. [37]. These authors analyzed monthly mean rainfall and observed that months with the lowest rainfall, from June to November, in the Sertão mesoregion in Paraíba State (also inserted in NEB), had the lowest variability. It is also worth noting the overestimation of estimated values in 2014 (N) and 2017 (MD) by the CHIRPS product. According to Rodrigues et al. [8], the CHIRPS product showed a moderate overestimation of monthly rainfall ranging between 0 and 100 mm while tending to underestimate values above 100 mm in the NEB.

3.3. Universal Kriging

Trend occurrence was detected using linear and quadratic polynomials, and for such cases, geostatistical analysis was carried out based on the residuals. Well-defined experimental semivariograms with high spatial dependence were produced based on the autocorrelation between the residual values for the years 1993, 1983, 2017, 2014, 2004, and 1985, as shown in Figure 5. From 1981 onwards, which is the initial date of availability of CHIRPS data, there was no SW year based on CHIRPS data (1981–2017). The semivariograms parameters and validation results are presented in Table 5. Figure 5 displays a set of modeled semivariograms for the six classes of SPI under study. The semivariances increase with the increase in the rainfall index, and therefore, the sill (C₀ + C₁) is directly related to the rainfall amount for the Brígida River basin, as also shown in Table 5.

The semivariogram models that best fit the data were the spherical (for ED and SD), followed by the Gaussian (for N, MW, and EW) and exponential (for MD) models. The relationship between the nugget effect and the sill ranged from 0.04% (MD) to 10.88% (EW), indicating a strong spatial dependence for all analyzed years, according to the classification of Cambardella et al. [59]. For monthly mean rainfall, [36] also found strong spatial dependence for most months.

Annual rainfall in the study area exhibited a spatial dependence structure, with ranges varying from 37.6 km (for EW) to 114.0 km (for MD), indicating a comprehensive structure of spatial dependence. Furthermore, Araújo et al. [27] also fitted semivariogram models for annual observed rainfall in the Brígida River basin, obtaining ranges from 20 to 45 km and a degree of spatial dependence ranging from moderate to strong. The coefficients of determination (R²) for the fitted models were moderate to high, with a minimum of 0.992 (for ED) and a maximum of 0.999 (for SD). Moreover, the semivariogram models were validated using the Jack-Knifing method. However, experimental semivariances presented gaps and oscillations due to uneven spatial distribution of rain gauges.

3.4. Cross-Semivariogram Models of Annual Rainfall

Well-constructed experimental cross-semivariograms were generated, exhibiting high spatial dependence based on the correlation between annual rainfall and elevation data for the years 1993, 1983, 2017, 2014, 2004, and 1985 (Figure 6). Parameters and validation results of the cross-semivariograms are shown in Table 6. The sill is generally related to the amount of rainfall, except for the year 1983 (SD), which can be explained by the lower value of the parameter C₀ (nugget effect).

3.5. Spatiotemporal Analysis (Rainfall UK, SGS and Cokriging)

The adjusted models were used to produce CHIRPS annual rainfall isoline maps for the years 1993 (ED), 1983 (SD), 2017 (MD), 2014 (N), 2004 (MW), and 1985 (EW) using universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK) techniques (Figure 7). The fitted models were also used to produce maps of annual rainfall with different runs (10, 100, and 1000) using the sequential Gaussian simulation (SGS) technique.

The northern region of the basin, which presents the highest elevation, encompasses the “Chapada do Araripe” highland, an area with considerable spatial variability of rainfall and a higher vegetation density. As expected, our study found that the amount of rainfall in the basin increased with elevation. Supporting the findings of this study, Mutiga et al. [67] calculated the water balance in the upper part of the northern Ewaso Ng’iro basin in Kenya, with altitudes ranging from 800 to 5200 m, and observed that with increasing elevation, the Normalized Difference Vegetation Index (NDVI) tended to increase due to greater rainfall occurrence, and actual evapotranspiration decreased. These results reinforce the well-known relationship between rainfall and elevation in hydrological systems.

The southern region of the Brígida River basin exhibited the lowest annual rainfall values for all years analyzed. This region is part of the Cabrobó desertification area located in the southern of the Pernambuco State and is characterized by a lack of vegetation cover and high degrees of soil erosion. Among the primary causes of the desertification process, low precipitation and deforestation should be noted [17]. The municipality of Exu (located at an altitude above 972 m, as shown in Figure 1b) consistently had the highest rainfall amounts in the basin for all maps.

The lowest rainfall occurred, in general, throughout the basin in the ED and SD years, while in the MD and N years, critical areas were located in the south and west. In the MW and EW years, the critical areas were in the southern part of the basin. The maps for the MD, N, MW, and EW years show the area of the Brígida basin that corresponds to the Köppen–Geiger climate classification of Aw (tropical climate with dry winter), as shown in Figure 1b.

Increasing the number of runs in the SGS maps often results in stability, with low or no improvement, and a convergence criterion can be achieved with 100 runs, reducing the computational effort, as observed in Bai and Tahmasebi [68] and Pereira et al. [69]. In addition, there was lower statistical error comparing the map of SGS-1000 with SGS-100. Furthermore, the SGS maps produced consistent results with low uncertainties and showed significant similarities to the maps of universal kriging (UK).

To evaluate the behavior of precipitation from the mouth to a high-altitude region in the basin, a transect analysis was performed on the maps of CHIRPS annual precipitation data produced using universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK), as shown in Figure 8.

3.6. Performance of Different Interpolation (UK, SGS, OCK) Methods for Annual Rainfall Estimation in the Study Area

The performance of the different interpolations (UK, SGS, OCK; Table 7) shows that SGS has the best performance, with a coefficient of determination (R²) of 100% and MBE and RMSE values close to zero. There is a slight advantage of OCK over UK with lower values of RMSE. This result is consistent with Adhikary et al. [39], who evaluated three geostatistical interpolation methods (ordinary kriging [OK], ordinary cokriging [OCK], external drift kriging [KED]) for mapping monthly rainfall from two basins in Australia and found that OCK interpolator yielded the best performance. In addition, maps of the standard deviations of universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK) are presented (Figure 9). It can be verified that uncertainties associated with the mapping are smaller when using the SGS, and the errors observed are more pronounced for the wetter years.

3.7. SPI Semivariogram Models

Rainfall data for each CHIRPS pixel was transformed into SPI data to investigate the spatial and temporal variation in historical rainfall patterns and drought severity. Since there is a spatial trend in the data, semivariograms were constructed with the residual values obtained for the years 1993, 1983, 2017, 2014, 2004, and 1985 (Figure 10), along with the adjustment parameters presented in Table 8.

The SPI showed a spatial dependence structure with ranges from 44.4 km (EW year) to 97.2 km (MD year), similar to the previous analysis with rainfall data. The semivariogram models that best fit the data were spherical (ED, SD, and EW years), Gaussian (N and MW), and exponential (MD). The relationship between the nugget effect and the sill of semivariograms ranged from 0.05% (EW) to 8.62% (N), indicating a strong spatial dependence for all analyzed years, according to the classification by Cambardella et al. [59]. The coefficients of determination (R²) for the model fitting were high, with a minimum of 0.992 (ED) and a maximum of 0.999 (N), similar to the results obtained previously with CHIRPS rainfall data. Furthermore, the semivariogram models were validated using the Jack-Knifing method.

3.8. SPI SGS

The fitted semivariogram models were used with the sequential Gaussian simulation (SGS) technique to produce maps of the SPI classifications with 100 runs (Figure 11). The municipalities of Ouricuri (ID-11), Parnamirim (ID-29), and Orocó (ID-34), as well as the northern regions of the municipalities of Santa Maria da Boa Vista (ID-35), Santa Cruz (ID-36), and Santa Filomena (ID-12), were identified as experiencing the most severe water deficits for all analyzed years, as classified by the SPI.

Even during Extremely Dry, Severely Dry, and Moderately Dry years, high (positive) SPI values are observed in the north of the hydrographic basin, which is a higher altitude area with an Aw (tropical climate with dry winter). This is consistent with a study conducted by Santos et al. [21], which estimated actual evapotranspiration (ETa) in the Brígida River basin and compared the dry year of 1990 and the wet year of 2009. The study found that areas with higher altitudes in the northern part of the basin had greater vegetation cover; lower surface temperature and albedo values; and higher NDVI (Normalized Difference Vegetation Index), SAVI (Soil-adjusted Vegetation Index), LAI (Leaf Area Index), and ETa (Actual evapotranspiration) values for both years analyzed. Therefore, a well-defined spatial pattern has been identified for these variables in the basin.

The anthropization of vegetated areas in the Caatinga biome is an extremely important component causing the reduction of precipitation, as reported by Silva et al. [70]. The authors found a reduction in precipitation for a deforested area compared to a region with native vegetation in the non-deforested Caatinga in the Pajeú watershed, Brazil. According to Panday et al. [71], the simultaneous occurrence of deforestation associated with climate variability has promoted significant impacts on the water balance in the Xingu watershed, Brazil, where abstraction of native vegetation significantly reduced rainfall in the basin. In the Brígida Basin, it was found that 57.32% of the area is anthropized and that the spatial distribution of rainfall is strongly influenced by the relief [20].

3.9. Analysis of Biophysical Indexes and Land Use and Land Cover (LULC)

Figure 12 shows the maps of biophysical indexes (Figure 12a,b) and major spatial classes of LULC (Figure 12c). In the maps for the MW years, it is possible to identify two reservoirs located in the central-west and central-south regions of the basin, which are the Entremontes Reservoir (capacity of 339.33 hm³, built in 1982) and Chapéu Reservoir (capacity of 188 hm³, built in 1936). These reservoirs are essential for domestic water supply, agricultural activities, and animal watering, according to the Brazilian National Water Authority [72,73]. Furthermore, the estimated area and stored volume of the reservoirs is consistent with the rainfall regime of a given year, assessed through the SPI-12 index (Table 9).

There is an increase in the area of water bodies and stored volume for the wetter rainfall regimes (from MD to EW). However, for the year 1985 (EW), there is a smaller area and volume than in 2004 (MW). This can be explained by the greater rainfall that occurred during the initial months of the wet season in 1985 (January and February), while in 2004, the highest rainfall occurred in the last month of the wet season (April). Additionally, this behavior is also associated with an increase in the use of stored water for irrigation. For a nearby basin, Silva et al. [24] estimated the total volume of water stored in small reservoirs to be 8.8 × 10⁵ m³, also using the MNDWI index. As expected, SPI-12 Index reflects long-term precipitation patterns and is linked to stored volumes in reservoirs due to the associated longer time scale.

The number of water bodies (blue pixels of the MNDWI) in the maps of the ED, MD, and N years is low (Figure 12a), as reservoir storage is strongly dependent of rainfall regime and water use allocation for different multiple uses. It should be highlighted that the hydrological river flow regime is typically temporary, and hence in periods of droughts, the river flow is interrupted until the following rainy season [74]. According to a study by Marengo et al. [28], which evaluated the climate patterns of drought in the NEB from 2010 to 2016, there was a sequence of dry years that affected the levels of reservoirs in the region and led to a water crisis. In addition, Silva et al. [23] evaluated the spatiotemporal dynamics of water bodies in the micro-region of Vale do Ipojuca (included in the NEB) from 2014 to 2021 and found a marked decrease in the area of water bodies between 2014 (1226 ha, the year of greatest water availability) and 2017 (683 ha), mainly in 2016 (558 ha). According to the authors, this dynamic is associated with the occurrence of severe droughts in the NEB, which caused significant reductions in water bodies volumes that lasted until mid-2019.

In relation to the total basin area (13,495.73 km²), 9.0%, 13.2%, 11.8%, and 31.2% were covered by arboreal Caatinga vegetation, considering a NDVI (Figure 12b) greater than or equal to 0.375, for the years MD (2017), N (2014), MW (2004), and EW (1985), respectively. These results support the findings of Marengo et al. [28], who evaluated the climatological aspects of the droughts that occurred in the semiarid region of northeast Brazil, and also concluded that the conditions of the La Niña event in 2013 were not enough to attenuate the drought conditions established years before, and the El Niño in 2015 worsened drought conditions in later years.

Similar results of NDVI were also obtained in studies by Brito et al. [13], who quantitatively evaluated drought events that occurred in the semiarid region of the NEB between 1981 and 2016 and found that droughts from 2011 to 2016 were more frequent, severe, and affected a larger area with significant impacts on the population and economic activities. Additionally, these results also support Cunha et al. [34], who reported that unlike other drought events that occurred in the NEB, the 2012–2013 drought persisted for 6 years (until 2017). One possible reason for the smaller areas corresponding to agriculture and livestock in the 2004 (MW) maps compared to the 2014 (N) maps may be due to the abandonment of agricultural areas from 2000 to 2009 in the semiarid region of Pernambuco State, influencing an expansion of exposed soil areas with low and stunted vegetation, as found by [22] while studying the changes caused by agricultural activities in the municipality of Capoeiras, also in the semiarid region of the Pernambuco State. A reduction in vegetation cover was also observed when comparing the years 1990 and 2009, according to a study by Santos et al. [21]. Furthermore, Silva et al. [23] evaluated the spatiotemporal dynamics of vegetation cover in the micro-region of Vale do Ipojuca (included in the NEB) from 2014 to 2021 and observed, based on the NDVI maps, that the year 2014 had denser vegetation and the resilience of the vegetation is observed in native vegetation (Caatinga) compared to 2017.

To the north of the basin, where Chapada do Araripe is located at higher elevations, there is a region with higher spatial variability of rainfall, as evidenced by the UK rainfall maps. This region also has greater vegetation density, as verified by the NDVI maps, which in turn influences the distribution of rainfall. It is noteworthy to compare Figure 11a–c (for SPI-12) to Figure 12b (MD year) in order to verify the rainfall–vegetation nexus in the basin, particularly in the higher areas. Even during dry years, the Caatinga canopy of arboreal species is present, highlighting the need to preserve and regulate land use in those areas. Undoubtedly, higher zones provide ecological services for the region, even during dry periods.

In their studies in the Brígida River basin, Galvíncio et al. [20] also found that the spatial distribution of rainfall is influenced by elevation. NDVI values between 0 and 0.250 are assumed to correspond to areas more susceptible to degradation, as they present areas with exposed soil and herbaceous-graminoid vegetation due to the low resilience of native vegetation in drier periods. Similar results were observed by Tomasella et al. [75], who monitored degraded areas in Northeast Brazil from 2000 to 2016 and found that NDVI values ranging from 0 to 0.30 were mostly associated with exposed soil and herbaceous-graminoid vegetation, as well as a very low density of subshrubs. Similar results were found by Silva et al. [22] when identifying this type of vegetation and exposed soil at NDVI values between 0.001 and 0.264. According to the authors, NDVI values lower than 0.264 were associated with areas occupied by agriculture or livestock.

In terms of the total area of the basin (13,495.73 km²), the Forest class (Caatinga vegetation) accounts for 66.6%, 66.6%, 69.6%, and 82.5% of it, according to LULC (Figure 12c), for the MD (2017), N (2014), MW (2004), and EW (1985) years, respectively. Meanwhile, the Farming class comprises 31.6%, 31.7%, 27.9%, and 15.7% of the basin for the same years. Indeed, there exists an inverse relationship between the Forest and Farming classes, where one class increases while the other decreases, mainly between the MD and EW years. For instance, the Forest class (Caatinga biome) had an area of 9155 and 11,335 km² for the MD and EW years, respectively, while the agricultural activity area accounted for 4339 and 2159 km² for the MD and EW years, respectively. During the transition period, there was a net increase of 2180 km² in the Forest class and a reduction of 2181 km² in the Farming class. This was due to the deforestation of Caatinga forest areas that were cleared for agricultural activities. These findings align with previous studies, where Silva et al. [22,23] evaluated degradation processes in the semiarid region of northeastern Brazil (NEB) using vegetation indices and LULC classes, highlighting the strong influence of agricultural activity on the Caatinga biome and vice versa.

The results of this study demonstrate the effectiveness of using spatial interpolation techniques and geostatistical analysis of rainfall data and provide valuable information for the management of public policies related to river basins at both national and international levels. These findings can aid in the characterization of river basins and inform the development of management and water management techniques.

4. Conclusions

We investigated rainfall variability patterns in a Brazilian semiarid basin and properly validated the CHIRPS rainfall product, which performed well in capturing rainfall and drought patterns for contrasting rainfall regimes. Among the semivariograms established and validated, the Gaussian model produced the best fit for most of the datasets that were used. Higher areas in the basin systematically present positive values for SPI-12, even in dry years, and exhibited high NDVI, which are directly related to canopy cover, ecological services, and streamflow regimes at the longer time scales. By using the sequential Gaussian simulation (SGS) technique, we identified low uncertainties in the CHIRPS data due to its proper reproduction of the spatial variability of annual rainfall. We suggest the adoption of 100 runs as a stopping criterion to reduce computational effort for simulations.

The use of sequential Gaussian simulation (SGS) for mapping both the annual precipitation and the Standard Precipitation Index (SPI-12) proved to be an efficient procedure, being able to provide essential data for hydrological analysis in semiarid regions. The developed framework successfully identified the vegetation cover dynamics and water features in the basin for different rainfall regimes. The Normalized Difference Vegetation Index (NDVI) properly detected the environmental degradation occurred during the 2012–2013 severe drought. Additionally, it is worth noting that the maps of generated SPI classifications provided valuable information on the spatial distribution of severe rainfall deficits in the basin. Through a spatiotemporal analysis of the indices, it was possible to identify the most vulnerable regions to intense droughts, thus emphasizing the need for effective management measures to be implemented for monitoring and mitigating their effects.

Moreover, the use of the Modified Normalized Difference Water Index (MNDWI) provided additional insights into water body dynamics in the basin and also indicated that, even in normal years, the number of water bodies is low, threatening water security in the region. Additionally, our study also analyzed the land use and land cover (LULC) in the basin, which revealed a significant change in the landscape over the years, with the expansion of cropland and livestock (Farming class) areas. The conversion of natural vegetation to Farming may have an impact on the water cycle and ecosystem services in the region, affecting both the water availability and quality. Effective management measures need to be implemented to maintain the water bodies and their associated ecological services. In conclusion, our study provides important insights into rainfall variability and hydrological processes in semiarid regions and presents a valuable framework for future studies in similar environments.

Author Contributions

Conceptualization, A.A.d.A.M., M.V.d.S. and L.d.B.d.S.; methodology, A.A.d.A.M., L.d.B.d.S., M.V.d.S., T.A.B.A. and A.A.d.C.; software, L.d.B.d.S., M.V.d.S., T.A.B.A. and A.A.d.C.; validation, L.d.B.d.S.; investigation, L.d.B.d.S. and M.V.d.S.; resources, L.d.B.d.S. and M.V.d.S.; data curation, L.d.B.d.S.; writing—original draft preparation, L.d.B.d.S.; writing—review and editing, A.A.d.A.M., J.L.M.P.d.L., M.V.d.S. and L.d.B.d.S.; visualization, A.A.d.A.M., T.G.F.d.S. and J.L.M.P.d.L.; supervision, A.A.d.A.M.; project administration, A.A.d.A.M. and L.d.B.d.S.; funding acquisition, A.A.d.A.M. and J.L.M.P.d.L. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data presented in this study are available upon request from the corresponding author.

Acknowledgments

We thank the Postgraduate Program in Agricultural Engineering (PGEA) of the Federal Rural University of Pernambuco (UFRPE) for supporting the development of this research.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Location of the São Francisco Basin, Brazil, South America (a); Brígida River basin, Pernambuco State (b); Shuttle Radar Topography Mission (SRTM) map (c); Köppen–Geiger climate classification and the spatial distribution of rainfall stations (ID) (d); and mean annual rainfall and years of record of the 40 analyzed rainfall stations (e).

Figure 2. Flowchart of the data analysis, processing steps, and procedures.

Figure 3. Correlation between observed annual rainfall at meteorological stations and simulated CHIRPS data (a) and cumulative density curve (b) for the period 1981–2017.

Figure 4. Box plots of CHIRPS annual rainfall data (C) and observed rainfall data (O) for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years.

Figure 5. Semivariogram models for CHIRPS residual annual rainfall data for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years. The separation distance is given in meters (m).

Figure 6. Cross-semivariogram models for CHIRPS annual rainfall and elevation data for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years. Separation distance is in meters (m).

Figure 7. Universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK) maps for CHIRPS annual rainfall data for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years, indicating the location of weather stations used in the study.

Figure 8. Transect analysis of CHIRPS annual precipitation data using universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK) maps for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years.

Figure 9. Standard deviation maps for universal kriging (UK), sequential Gaussian simulation (SGS), and ordinary cokriging (OCK) of CHIRPS annual rainfall data. In addition, errors (SGS-UK) and (SGS-OCK) for Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004) and Extremely Wet (1985) years.

Figure 10. Semivariogram models fitted for residual SPI data during Extremely Dry (1993), Severely Dry (1983), Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years. Separation distance is in meters (m).

Figure 11. Sequential Gaussian simulation (SGS) maps for the SPI with 100 runs for Extremely Dry (1993) (a), Severely Dry (1983) (b), Moderately Dry (2017) (c), Normal (2014) (d), Moderately Wet (2004) (e), and Extremely Wet (1985) (f) years, indicating the location of weather stations used in the study.

Figure 12. Modified Normalized Difference Water Index (MNDWI) maps (a), Normalized Difference Vegetation Index (NDVI) maps (b), and annual spatiotemporal dynamics of land use and land cover (LULC) maps (c) for Moderately Dry (2017), Normal (2014), Moderately Wet (2004), and Extremely Wet (1985) years, indicating the location of the weather stations used in the study.

Table 1

SPI classification according to McKee et al. [11].

SPI	Classifications
≥2.00	Extremely Wet (EW)
1.00 to 1.99	Severely Wet (SW)
0.50 to 0.99	Moderately Wet (MW)
0.49 to −0.49	Normal (N)
−0.50 to −0.99	Moderately Dry (MD)
−1.00 to −1.99	Severely Dry (SD)
≤−2.00	Extremely Dry (ED)

Table 2

Satellite images used in the study with hour and date.

Year	Satellite	Date	Hour (UTM)	Orbit	Point
ED (1993)	No Landsat images after September and no clouds
SD (1983)	No Landsat images
MD (2017)	Landsat 8	30 September 2017	12:47:34	217	65
		30 September 2017	12:47:58	217	66
		7 October 2017	12:53:47	218	65
N (2014)	Landsat 8	22 September 2014	12:47:26	217	65
		22 September 2014	12:47:50	217	66
		13 September 2014	12:53:40	218	65
MW (2004)	Landsat 5	13 November 2004	12:32:25	217	65
		13 November 2004	12:32:49	217	66
		4 November 2004	12:38:28	218	65
EW (1985)	Landsat 5	24 October 1985	12:16:14	217	65
		24 October 1985	12:16:38	217	66
		15 October 1985	12:22:31	218	65

Source: [63].

Table 3

SPI classification for the total annual rainfall measured at the station locations.

Year	SPI	Classification	Year	SPI	Classification	Year	SPI	Classification	Year	SPI	Classification
1962	0.47	N	1976	−0.33	N	1990	−1.27	MD	2004	1.32	MW
1963	0.55	N	1977	0.42	N	1991	−0.81	N	2005	0.29	N
1964	1.72	SW	1978	−0.04	N	1992	−0.78	N	2006	−0.32	N
1965	0.19	N	1979	0.70	N	1993	−2.21	ED	2007	−0.62	N
1966	0.11	N	1980	0.49	N	1994	0.01	N	2008	0.59	N
1967	1.09	MW	1981	−0.19	N	1995	0.29	N	2009	0.94	N
1968	0.29	N	1982	−1.01	MD	1996	0.68	N	2010	−0.15	N
1969	−0.13	N	1983	−1.98	SD	1997	0.54	N	2011	0.41	N
1970	−0.48	N	1984	0.70	N	1998	−1.26	MD	2012	−2.47	ED
1971	0.50	N	1985	3.30	EW	1999	0.05	N	2013	−0.78	N
1972	−0.65	N	1986	0.35	N	2000	0.02	N	2014	−0.27	N
1973	0.67	N	1987	−0.80	N	2001	−0.86	N	2015	−0.37	N
1974	1.84	SW	1988	1.31	MW	2002	−0.37	N	2016	−1.32	MD
1975	0.12	N	1989	1.25	MW	2003	−0.70	N	2017	−1.02	MD

Where ED = Extremely Dry, SD = Severely Dry, MD = Moderately Dry, N = Normal, MW = Moderately Wet, SW = Severely Wet, EW = Extremely Wet.

Table 4

Descriptive statistics of annual rainfall estimated by CHIRPS for the Brígida basin.

Year	Mean	Median	Maximum	Minimum	n	SD	CV (%)	p-Value (KS)
Extremely Dry (1993)	329.48	280.82	843.10	138.21	930	137.17	41.63%	0.14
Severely Dry (1983)	361.05	359.68	658.09	175.30	930	76.19	21.10%	0.05
Moderately Dry (2017)	475.55	445.60	1001.36	212.67	930	145.69	30.64%	0.10
Normal (2014)	641.66	543.09	1444.71	333.15	930	248.55	38.74%	0.16
Moderately Wet (2004)	854.19	824.72	1566.14	355.03	930	213.5	24.99%	0.09
Extremely Wet (1985)	1292.89	1237.86	2103.26	692.42	930	243.65	18.88%	0.09

Where SD = standard deviation, CV = coefficient of variation, n = number of sampling points.

Table 5

Semivariogram parameters and models for annual rainfall residual data; statistics of Jack-Knifing procedure.

SPI Class	Model	Nugget (C0)	Sill (C0 + C1)	Range (A, km)	R2	DSD (%)	Jack-Knifing
							M	SD
Extremely Dry (1993)	Spherical	80	1348	55.9	0.992	5.93%	0.01	1.04
Severely Dry (1983)	Spherical	59	1998	61.0	0.999	2.95%	0.01	1.03
Moderately Dry (2017)	Exponential	1	2405	114.0	0.998	0.04%	0.01	1.02
Normal (2014)	Gaussian	580	6059	50.2	0.998	9.57%	0.01	1.06
Moderately Wet (2004)	Gaussian	370	7852	55.4	0.996	4.71%	0.01	1.04
Extremely Wet (1985)	Gaussian	1160	10,660	37.6	0.996	10.88%	0.01	1.04

Where DSD = degree of spatial dependence, M = mean, SD = standard deviation.

Table 6

Cross-semivariogram parameters and models for annual rainfall and elevation; statistics of Jack-Knifing procedure.

SPI Class	Model	Nugget (C0)	Sill (C0 + C1)	Range (A, km)	R2	DSD (%)	Jack-Knifing
							M	SD
Extremely Dry (1993)	Gaussian	100	32,040	168.5	0.987	0.31%	0.00	1.00
Severely Dry (1983)	Gaussian	10	17,200	232.3	0.986	0.06%	0.01	1.01
Moderately Dry (2017)	Gaussian	100	37,160	217.0	0.985	0.27%	0.00	1.00
Normal (2014)	Gaussian	100	44,920	167.5	0.979	0.22%	0.00	1.00
Moderately Wet (2004)	Gaussian	100	51,300	212.9	0.990	0.19%	0.00	1.00
Extremely Wet (1985)	Gaussian	100	50,690	171.9	0.979	0.20%	0.00	1.00

Where DSD = degree of spatial dependence, M = mean, SD = standard deviation.

Table 7

Performance of different interpolation methods (UK, SGS, OCK) for annual rainfall estimation in the study area.

SPI Class	MBE (mm)			RMSE (mm)			R²
	UK	SGS	OCK	UK	SGS	OCK	UK	SGS	OCK
Extremely Dry (1993)	−0.1415	−0.0003	−0.1615	55.877	0.0029	30.632	0.9972	1.0000	0.9992
Severely Dry (1983)	0.0010	−0.0001	0.0205	57.187	0.0028	36.781	0.9895	1.0000	0.9955
Moderately Dry (2017)	0.6475	0.0000	0.6481	35.520	0.0029	35.248	0.9984	1.0000	0.9984
Normal (2014)	0.1050	−0.0001	0.3873	18.392	0.0029	54.153	0.9910	1.0000	0.9992
Moderately Wet (2004)	0.3519	0.0000	−0.0886	160.142	0.0029	47.650	0.9896	1.0000	0.9991
Extremely Wet (1985)	0.4379	−0.0003	0.7608	309.794	0.0028	87.278	0.9770	1.0000	0.9982

Where MBE = mean bias error, RMSE = root mean squared error, and R² = coefficient of determination.

Table 8

Semivariogram adjustment parameters of SPI residual data; Jack-Knifing statistics.

SPI Class	Model	Nugget (C0)	Sill (C0 + C1)	Range (A, km)	R2	DSD (%)	Jack-Knifing
							M	SD
Extremely Dry (1993)	Spherical	0.0044	0.0833	54.2	0.992	5.28%	0.01	1.04
Severely Dry (1983)	Spherical	0.0180	0.3590	60.4	0.998	5.01%	0.01	1.04
Moderately Dry (2017)	Exponential	0.0001	0.1078	97.2	0.995	0.09%	0.01	1.03
Normal (2014)	Gaussian	0.0099	0.1148	47.8	0.999	8.62%	0.01	1.06
Moderately Wet (2004)	Gaussian	0.0110	0.1710	53.7	0.997	6.43%	0.01	1.05
Extremely Wet (1985)	Spherical	0.0001	0.1832	44.4	0.995	0.05%	0.01	1.03

Where DSD = degree of spatial dependence, M = mean, SD = standard deviation.

Table 9

Water mirror area and volume of water stored in the water reservoirs of the Brígida River basin, semiarid region of Pernambuco, Brazil.

Year	Total Area of Entremontes Reservoir (m² × 105)	Total Area of Chapéu Reservoir (m² × 105)	Total Volume of Water Stored in EntremontesReservoir (m³ × 105)	Total Volume of Water Stored in ChapéuReservoir (m³ × 105)
Moderately Dry (2017)	3.35	0	6.71	0
Normal (2014)	38.0	8.28	76.4	16.6
Moderately Wet (2004)	327	191	659	386
Extremely Wet (1985)	259	1.47	522	2.94

There are no images for 1983 and 1993.

References

1. Cristiano, E.; ten Veldhuis, M.C.; van de Giesen, N. Spatial and temporal variability of rainfall and their effects on hydrological response in urban areas—A review. Hydrol. Earth Syst. Sci.; 2017; 21, pp. 3859-3878. [DOI: https://dx.doi.org/10.5194/hess-21-3859-2017]

2. Peña-Angulo, D.; Nadal-Romero, E.; González-Hidalgo, J.C.; Albaladejo, J.; Andreu, V.; Bagarello, V.; Barhi, H.; Batalla, R.J.; Bernal, S.; Bienes, R. et al. Spatial variability of the relationships of runoff and sediment yield with weather types throughout the Mediterranean basin. J. Hydrol.; 2019; 571, pp. 390-405. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2019.01.059]

3. Kim, J.; Lee, J.; Kim, D.; Kang, B. The role of rainfall spatial variability in estimating areal reduction factors. J. Hydrol.; 2019; 568, pp. 416-426. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2018.11.014]

4. Peña-Angulo, D.; Nadal-Romero, E.; González-Hidalgo, J.C.; Albaladejo, J.; Andreu, V.; Barhi, H.; Bernal, S.; Biddoccu, M.; Bienes, R.; Campo, J. et al. Relationship of weather types on the seasonal and spatial variability of rainfall, runoff, and sediment yield in the western Mediterranean basin. Atmosphere; 2020; 11, 609. [DOI: https://dx.doi.org/10.3390/atmos11060609]

5. Da Silva, J.L.B.; Moura, G.B.D.A.; da Silva, M.V.; Lopes, P.M.O.; Guedes, R.V.D.S.; Silva, F.D.F.E.; Ortiz, P.F.S.; Rodrigues, J.A.D.M. Changes in the water resources, soil use and spatial dynamics of Caatinga vegetation cover over semiarid region of the Brazilian Northeast. Remote Sens. Appl. Soc. Environ.; 2020; 20, 100372. [DOI: https://dx.doi.org/10.1016/j.rsase.2020.100372]

6. Santos, J.F.S.; Naval, L.P. Spatial and temporal dynamics of water footprint for soybean production in areas of recent agricultural expansion of the Brazilian savannah (Cerrado). J. Clean. Prod.; 2020; 251, 119482. [DOI: https://dx.doi.org/10.1016/j.jclepro.2019.119482]

7. Singh, V.P. Handbook of Applied Hydrology; 2nd ed. McGraw-Hill Education—Europe: London, UK, 2016; ISBN 0071835091

8. Paredes-Trejo, F.J.; Barbosa, H.A.; Kumar, T.V.L. Validating CHIRPS-based satellite precipitation estimates in Northeast Brazil. J. Arid Environ.; 2017; 139, pp. 26-40. [DOI: https://dx.doi.org/10.1016/j.jaridenv.2016.12.009]

9. Rodrigues, D.T.; Gonçalves, W.A.; Spyrides, M.H.C.; e Silva, C.M.S. Spatial and temporal assessment of the extreme and daily precipitation of the Tropical Rainfall Measuring Mission satellite in Northeast Brazil. Int. J. Remote Sens.; 2019; 41, pp. 549-572. [DOI: https://dx.doi.org/10.1080/01431161.2019.1643940]

10. Neto, R.M.B.; Santos, C.A.G.; Silva, J.F.C.B.D.C.; da Silva, R.M.; dos Santos, C.A.C.; Mishra, M. Evaluation of the TRMM product for monitoring drought over Paraíba State, northeastern Brazil: A trend analysis. Sci. Rep.; 2021; 11, 1097. [DOI: https://dx.doi.org/10.1038/s41598-020-80026-5]

11. McKee, T.B.; Doesken, N.J.; Kleist, J. The Relationship of Drought Frequency and Duration to Time Scales. Proceedings of the 8th Conference on Applied Climatology; Anaheim, CA, USA, 17–22 January 1993; pp. 179-184.

12. Juliani, B.H.T.; Okawa, C.M.P. Application of a standardized precipitation index for meteorological drought analysis of the semi-arid climate influence in Minas Gerais, Brazil. Hydrology; 2017; 4, 26. [DOI: https://dx.doi.org/10.3390/hydrology4020026]

13. Brito, S.S.B.; Cunha, A.P.M.A.; Cunningham, C.C.; Alvalá, R.C.; Marengo, J.A.; Carvalho, M.A. Frequency, duration and severity of drought in the Semiarid Northeast Brazil region. Int. J. Climatol.; 2018; 38, pp. 517-529. [DOI: https://dx.doi.org/10.1002/joc.5225]

14. Santos, C.A.G.; Brasil Neto, R.M.; de Araújo Passos, J.S.; da Silva, R.M. Drought assessment using a TRMM-derived standardized precipitation index for the upper São Francisco River basin, Brazil. Environ. Monit. Assess.; 2017; 189, 250. [DOI: https://dx.doi.org/10.1007/s10661-017-5948-9]

15. Tsesmelis, D.E.; Leveidioti, I.; Karavitis, C.A.; Kalogeropoulos, K.; Vasilakou, C.G.; Tsatsaris, A.; Zervas, E. Spatiotemporal Application of the Standardized Precipitation Index (SPI) in the Eastern Mediterranean. Climate; 2023; 11, 95. [DOI: https://dx.doi.org/10.3390/cli11050095]

16. Toté, C.; Patricio, D.; Boogaard, H.; van der Wijngaart, R.; Tarnavsky, E.; Funk, C. Evaluation of satellite rainfall estimates for drought and flood monitoring in Mozambique. Remote Sens.; 2015; 7, pp. 1758-1776. [DOI: https://dx.doi.org/10.3390/rs70201758]

17. Mariano, D.A.; dos Santos, C.A.; Wardlow, B.D.; Anderson, M.C.; Schiltmeyer, A.V.; Tadesse, T.; Svoboda, M.D. Use of remote sensing indicators to assess effects of drought and human-induced land degradation on ecosystem health in Northeastern Brazil. Remote Sens. Environ.; 2018; 213, pp. 129-143. [DOI: https://dx.doi.org/10.1016/j.rse.2018.04.048]

18. Paca, V.H.d.M.; Espinoza-Dávalos, G.E.; Moreira, D.M.; Comair, G. Variability of Trends in Precipitation across the Amazon River Basin Determined from the CHIRPS Precipitation Product and from Station Records. Water; 2020; 12, 1244. [DOI: https://dx.doi.org/10.3390/w12051244]

19. Rouse, J.W.; Haas, R.H.; Schell, J.A.; Deering, D.W.; Harlan, J.C. Monitoring the Vernal Advancement of Retrogradation of Natural Vegetation; NASA: Washington, DC, USA, 1973.

20. Galvíncio, J.D.; Sá, I.I.S.; de Moura, M.S.B.; Ribeiro, J.G. Determinação das características físicas, climáticas e da paisagem da bacia hidrográfica do rio Brigida com o auxilio de técnicas de geoprocessamento e sensoriamento remoto. Rev. Geogr.; 2007; 24, pp. 83-96.

21. Santos, C.A.G.; da Silva, R.M.; Silva, A.M.; Neto, R.M.B. Estimation of evapotranspiration for different land covers in a Brazilian semi-arid region: A case study of the Brígida River basin, Brazil. J. S. Am. Earth Sci.; 2017; 74, pp. 54-66. [DOI: https://dx.doi.org/10.1016/j.jsames.2017.01.002]

22. Da Silva, M.V.; Pandorfi, H.; Lopes, P.M.O.; da Silva, J.L.B.; de Almeida, G.L.P.; Silva, D.A.D.O.; dos Santos, A.; Rodrigues, J.A.D.M.; Batista, P.H.D.; Jardim, A.M.D.R.F. Pilot monitoring of caatinga spatial-temporal dynamics through the action of agriculture and livestock in the brazilian semiarid. Remote Sens. Appl. Soc. Environ.; 2020; 19, 100353. [DOI: https://dx.doi.org/10.1016/j.rsase.2020.100353]

23. Da Silva, M.V.; Pandorfi, H.; de Oliveira-Júnior, J.F.; da Silva, J.L.B.; de Almeida, G.L.P.; Montenegro, A.A.d.A.; Mesquita, M.; Ferreira, M.B.; Santana, T.C.; Marinho, G.T.B. et al. Remote sensing techniques via Google Earth Engine for land degradation assessment in the Brazilian semiarid region, Brazil. J. S. Am. Earth Sci.; 2022; 120, 104061. [DOI: https://dx.doi.org/10.1016/j.jsames.2022.104061]

24. Silva, J.R.I.; Montenegro, A.A.D.A.; Farias, C.W.L.D.A.; Jardim, A.M.D.R.F.; da Silva, T.G.F.; Montenegro, S.M.G.L. Morphometric characterization and land use of the Pajeú river basin in the Brazilian semi-arid region. J. S. Am. Earth Sci.; 2022; 118, 103939. [DOI: https://dx.doi.org/10.1016/j.jsames.2022.103939]

25. De Lima, I.P.; Jorge, R.G.; de Lima, J.L.M.P. Remote Sensing Monitoring of Rice Fields: Towards Assessing Water Saving Irrigation Management Practices. Front. Remote Sens.; 2021; 2, 762093. [DOI: https://dx.doi.org/10.3389/frsen.2021.762093]

26. Silva, A.P.; de Lima, I.P.; Santo, M.F.E.S.; Pires, V.C. Assessing changes in drought and wetness episodes in drainage basins using the Standardized Precipitation Index. Die Bodenkult.; 2014; 65, pp. 31-37.

27. Araújo, H.L.; Montenegro, A.A.d.A.; Lopes, I.; de Carvalho, A.A.; e Silva, E.C.; Gonçalves, G.E. Precipitation spatialization in the Brígida river basin in Pernambuco semi-arid. Rev. Bras. Geogr. Fis.; 2020; 13, pp. 391-405. [DOI: https://dx.doi.org/10.26848/rbgf.v13.1.p391-405]

28. Marengo, J.A.; Alves, L.M.; Alvala, R.C.S.; Cunha, A.P.; Brito, S.; Moraes, O.L.L. Climatic characteristics of the 2010-2016 drought in the semiarid northeast Brazil region. An. Acad. Bras. Cienc.; 2018; 90, pp. 1973-1985. [DOI: https://dx.doi.org/10.1590/0001-3765201720170206]

29. Montenegro, A.A.A.; Ragab, R. Hydrological response of a Brazilian semi-arid catchment to different land use and climate change scenarios: A modelling study. Hydrol. Process.; 2010; 24, pp. 2705-2723. [DOI: https://dx.doi.org/10.1002/hyp.7825]

30. Wanderley, H.S.; de Amorim, R.F.C.; de Carvalho, F.O. Variabilidade espacial e preenchimento de falhas de dados pluviométricos para o estado de Alagoas. Rev. Bras. Meteorol.; 2012; 27, pp. 347-354. [DOI: https://dx.doi.org/10.1590/S0102-77862012000300009]

31. Hengl, T.; Walsh, M.G.; Sanderman, J.; Wheeler, I.; Harrison, S.P.; Prentice, I.C. Global mapping of potential natural vegetation: An assessment of machine learning algorithms for estimating land potential. PeerJ; 2018; 2018, e5457. [DOI: https://dx.doi.org/10.7717/peerj.5457]

32. Wolff, W.; Duarte, S.N.; Fernandes, R.D.M. Spatialization of the annual and seasonal average precipitations in the state of Santa Catarina, Brazil. J. S. Am. Earth Sci.; 2020; 103, 102735. [DOI: https://dx.doi.org/10.1016/j.jsames.2020.102735]

33. Da Silva, M.V.; Pandorfi, H.; Jardim, A.M.D.R.F.; de Oliveira-Júnior, J.F.; da Divincula, J.S.; Giongo, P.R.; da Silva, T.G.F.; de Almeida, G.L.P.; Moura, G.B.D.A.; Lopes, P.M.O. Spatial modeling of rainfall patterns and groundwater on the coast of northeastern Brazil. Urban Clim.; 2021; 38, 100911. [DOI: https://dx.doi.org/10.1016/j.uclim.2021.100911]

34. Cunha, A.P.M.A.; Tomasella, J.; Ribeiro-Neto, G.G.; Brown, M.; Garcia, S.R.; Brito, S.B.; Carvalho, M.A. Changes in the spatial–temporal patterns of droughts in the Brazilian Northeast. Atmos. Sci. Lett.; 2018; 19, e855. [DOI: https://dx.doi.org/10.1002/asl.855]

35. Da Silva, R.M.; Santos, C.A.G.; Silva, J.F.C.B.D.C.; Silva, A.M.; Neto, R.M.B. Spatial distribution and estimation of rainfall trends and erosivity in the Epitácio Pessoa reservoir catchment, Paraíba, Brazil. Nat. Hazards; 2020; 102, pp. 829-849. [DOI: https://dx.doi.org/10.1007/s11069-020-03926-9]

36. Medauar, C.C.; Silva, S.A.; Carvalho, L.C.C.; Galvão, Í.M.; Macêdo, P.V. Spatial-temporal variability of rainfall and mean air temperature for the state of Bahia, Brazil. An. Acad. Bras. Cienc.; 2020; 92, e20181283. [DOI: https://dx.doi.org/10.1590/0001-3765202020181283]

37. De Medeiros, E.S.; de Lima, R.R.; de Olinda, R.A.; dos Santos, C.A.C. Modeling spatiotemporal rainfall variability in Paraíba, Brazil. Water; 2019; 11, 1843. [DOI: https://dx.doi.org/10.3390/w11091843]

38. Page, T.; Beven, K.J.; Hankin, B.; Chappell, N.A. Interpolation of rainfall observations during extreme rainfall events in complex mountainous terrain. Hydrol. Process.; 2022; 36, 14758. [DOI: https://dx.doi.org/10.1002/hyp.14758]

39. Adhikary, S.K.; Muttil, N.; Yilmaz, A.G. Cokriging for enhanced spatial interpolation of rainfall in two Australian catchments. Hydrol. Process.; 2017; 31, pp. 2143-2161. [DOI: https://dx.doi.org/10.1002/hyp.11163]

40. Deutsch, C.V.; Journel, A.G. GSLIB: Geostatistical Software Library and User’s Guide; 2nd ed. Oxford University Press: New York, NY, USA, 1998; ISBN 0195073924

41. Méndez-Venegas, J.; Díaz-Viera, M.A.; Herrera, G.S.; Valdés-Manzanilla, A. Geostatistical simulation of spatial variability of convective storms in Mexico City Valley. Geofis. Int.; 2013; 52, pp. 111-120. [DOI: https://dx.doi.org/10.1016/S0016-7169(13)71466-1]

42. Wang, G.; Gertner, G.; Singh, V.; Shinkareva, S.; Parysow, P.; Anderson, A. Spatial and Temporal Prediction and Uncertainty Analysis of Rainfall Erosivity for the Revised Universal Soil Loss Equation. Ecol. Modell.; 2002; 153, pp. 143-155. [DOI: https://dx.doi.org/10.1016/S0304-3800(01)00507-5]

43. Goovaerts, P. Geostatistics for Natural Resources Evaluation; Oxford University Press: New York, NY, USA, 1997; ISBN 0195115384

44. Isaaks, E.H.; Srivastava, R.M. An Introduction to Applied Geostatistics; Oxford University Press: New York, NY, USA, 1989; ISBN 0195050134

45. Goovaerts, P. Geostatistical approaches for incorporating elevation into the spatial interpolation of rainfall. J. Hydrol.; 2000; 228, pp. 113-129. [DOI: https://dx.doi.org/10.1016/S0022-1694(00)00144-X]

46. Cuartas, L.A.; Cunha, A.P.M.D.A.; Alves, J.A.; Parra, L.M.P.; Deusdará-Leal, K.; Costa, L.C.O.; Molina, R.D.; Amore, D.; Broedel, E.; Seluchi, M.E. et al. Recent Hydrological Droughts in Brazil and Their Impact on Hydropower Generation. Water; 2022; 14, 601. [DOI: https://dx.doi.org/10.3390/w14040601]

47. Agência Pernambucana de Águas e Clima (APAC). Relatório de Situação de Recursos Hídricos do Estado de Pernambuco 2011/2012. Available online: https://www.lai.pe.gov.br/apac/wp-content (accessed on 16 June 2022).

48. Instituto de Pesquisa Econômica Aplicada (IPEA). Transposição do Rio São Francisco: Análise de Oportunidade do Projeto; IPEA: Rio de Janeiro, Brazil, 2011; ISBN 1415-4765

49. Alvares, C.A.; Stape, J.L.; Sentelhas, P.C.; Gonçalves, J.L.D.M.; Sparovek, G. Köppen’s climate classification map for Brazil. Meteorol. Z.; 2013; 22, pp. 711-728. [DOI: https://dx.doi.org/10.1127/0941-2948/2013/0507] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/24622815]

50. Beck, H.E.; Zimmermann, N.E.; McVicar, T.R.; Vergopolan, N.; Berg, A.; Wood, E.F. Present and future köppen-geiger climate classification maps at 1-km resolution. Sci. Data; 2018; 5, 180214. [DOI: https://dx.doi.org/10.1038/sdata.2018.214] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30375988]

51. Instituto Nacional de Meteorologia (INMET). Normais Climatológicas. Available online: https://clima.inmet.gov.br/GraficosClimatologicos/ (accessed on 10 June 2022).

52. Soil Survey Staff. Keys to Soil Taxonomy; 9th ed. US Department of Agriculture: Washington, DC, USA, 2006.

53. Gerardo, R.; de Lima, I.P. Sentinel-2 Satellite Imagery-Based Assessment of Soil Salinity in Irrigated Rice Fields in Portugal. Agriculture; 2022; 12, 1490. [DOI: https://dx.doi.org/10.3390/agriculture12091490]

54. Agência Pernambucana de Águas e Clima (APAC). Coordenadas das estações pluviométricas do Estado de Pernambuco. Available online: http://old.apac.pe.gov.br/meteorologia/coordenadas/index.html (accessed on 16 June 2022).

55. Warrick, A.W.; Nielsen, D.R. Spatial variability of soil physical properties in the field. Applications of Soil Physics; Academic Press: New York, NY, USA, 1980.

56. Journel, A.G. Fundamentals of Geostatistics in Five Lessons; American Geophysical Union: Washington, DC, USA, 1989; ISBN 0-87590-708-3

57. Landim, P.M.B. Análise Estatística de Dados Geológicos; 2nd ed. Edunesp: Rio Claro, Brazil, 2003; ISBN 9788571395046

58. Journel, A.G.; Huijbregts, C.J. Mining Geostatistics; Oxford University Press: London, UK, 1978; ISBN 0123910501

59. Cambardella, C.A.; Moorman, T.B.; Novak, J.M.; Parkin, T.B.; Karlen, D.L.; Turco, R.F.; Konopka, A.E. Field-Scale Variability of Soil Properties in Central Iowa Soils. Soil Sci. Soc. Am. J.; 1994; 58, pp. 1501-1511. [DOI: https://dx.doi.org/10.2136/sssaj1994.03615995005800050033x]

60. Vauclin, M.; Vieira, S.R.; Vachaud, G.; Nielsen, D.R. The Use of Cokriging with Limited Field Soil Observations. Soil Sci. Soc. Am. J.; 1983; 47, pp. 175-184. [DOI: https://dx.doi.org/10.2136/sssaj1983.03615995004700020001x]

61. Hengl, T.; Heuvelink, G.B.M.; Stein, A. Comparison of Kriging with External Drift and Regression-Kriging; ITC: Geneve, Switzerland, 2003; 17.

62. Goovaerts, P. Geostatistical modelling of uncertainty in soil science. Geoderma; 2001; 103, pp. 3-26. [DOI: https://dx.doi.org/10.1016/S0016-7061(01)00067-2]

63. Google Earth Engine. Landsat Collections. Available online: https://developers.google.com/earth-engine/datasets/catalog/landsat (accessed on 10 May 2022).

64. Xu, H. Modification of normalised difference water index (NDWI) to enhance open water features in remotely sensed imagery. Int. J. Remote Sens.; 2006; 27, pp. 3025-3033. [DOI: https://dx.doi.org/10.1080/01431160600589179]

65. Molle, F. Geometria de Pequenos Açudes; SUDENE: Recife, Brazil, 1994; Volume 1, ISBN 55655812813

66. MapBiomas Brazil. Coleção 7.1 da Série Anual de Mapas de Cobertura e Uso de Solo do Brasil. Available online: https://mapbiomas.org/colecoes-mapbiomas-1?cama_set_language=pt-BR (accessed on 11 March 2023).

67. Mutiga, J.K.; Su, Z.; Woldai, T. Using satellite remote sensing to assess evapotranspiration: Case study of the upper Ewaso Ng’iro North Basin, Kenya. Int. J. Appl. Earth Obs. Geoinf.; 2010; 12, pp. 100-108. [DOI: https://dx.doi.org/10.1016/j.jag.2009.09.012]

68. Bai, T.; Tahmasebi, P. Sequential Gaussian simulation for geosystems modeling: A machine learning approach. Geosci. Front.; 2022; 13, 101258. [DOI: https://dx.doi.org/10.1016/j.gsf.2021.101258]

69. Pereira, G.T.; Teixeira, D.D.B.; de Souza, Z.M.; de Oliveira, I.R.; Marques Júnior, J. Stochastic simulations of calcium contents in sugarcane area. Rev. Bras. Eng. Agric. Ambient.; 2015; 19, pp. 767-772. [DOI: https://dx.doi.org/10.1590/1807-1929/agriambi.v19n8p767-772]

70. Da Silva, T.G.F.; de Queiroz, M.G.; Zolnier, S.; de Souza, L.S.B.; de Souza, C.A.A.; de Moura, M.S.B.; de Araújo, G.G.L.; Neto, A.J.S.; dos Santos, T.S.; de Melo, A.L. et al. Soil properties and microclimate of two predominant landscapes in the Brazilian semiarid region: Comparison between a seasonally dry tropical forest and a deforested area. Soil Tillage Res.; 2021; 207, 104852. [DOI: https://dx.doi.org/10.1016/j.still.2020.104852]

71. Panday, P.K.; Coe, M.T.; Macedo, M.N.; Lefebvre, P.; de Castanho, A.D. Deforestation offsets water balance changes due to climate variability in the Xingu River in eastern Amazonia. J. Hydrol.; 2015; 523, pp. 822-829. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2015.02.018]

72. Agência Nacional de Águas e Saneamento Básico (ANA). Reservatórios do Semiárido Brasileiro: Hidrologia, Balanço Hídrico e Operação—Reservatório Entremontes. Available online: http://portal1.snirh.gov.br/arquivos/semiarido/204res/SF_Entremontes.pdf (accessed on 1 February 2022).

73. Agência Nacional de Águas e Saneamento Básico (ANA). Reservatórios do Semiárido Brasileiro: Hidrologia, Balanço Hídrico e Operação—Reservatório Chapéu. Available online: http://portal1.snirh.gov.br/arquivos/semiarido/204res/SF_Chapéu.pdf (accessed on 1 February 2022).

74. De Araújo, J.C.; Mamede, G.L.; De Lima, B.P. Hydrological guidelines for reservoir operation to enhance water governance: Application to the Brazilian Semiarid region. Water; 2018; 10, 1628. [DOI: https://dx.doi.org/10.3390/w10111628]

75. Tomasella, J.; Silva Pinto Vieira, R.M.; Barbosa, A.A.; Rodriguez, D.A.; de Oliveira Santana, M.; Sestini, M.F. Desertification trends in the Northeast of Brazil over the period 2000–2016. Int. J. Appl. Earth Obs. Geoinf.; 2018; 73, pp. 197-206. [DOI: https://dx.doi.org/10.1016/j.jag.2018.06.012]