DOI: 10.20937/ATM.53103
Received: July 15, 2021; Accepted: April 25, 2022
Javier de Jesus CorreaIslas^{1}, Juan Manuel RomeroPadilla^{1}*, Paulino PérezRodríguez^{1 }and Antonio VázquezAlarcón^{2}
^{1} PSEIEstadística, Colegio de Postgraduados Campus Montecillo, Km. 36.5, MéxicoTexcoco, 56230 Montecillo, Estado de México.
^{2} Departamento de Suelos, Universidad Autónoma Chapingo, km 38.5 Carretera MéxicoTexcoco, Chapingo, 56230, Estado de México.
* Corresponding author: romero.manuel@colpos.mx
RESUMEN
Se calculó un mapa de temperatura media anual mediante el método de interpolación de Kriging para la zona centronorte de México con la finalidad de obtener la condición actual de aridez, así como posibles escenarios en el futuro cercano y lejano. El gradiente altitudinal se estimó mediante regresión lineal y se usó en la estimación de la temperatura media. Las Áreas de Influencia Climática (CIA) se obtuvieron superponiendo la capa de precipitación oficial y la capa de temperatura media anual con la ayuda de herramientas de Sistemas de Información Geográfica. Se generaron bases de datos mensuales de variables climáticas para cada CIA y se estimó la evapotranspiración potencial utilizando la metodología de Thorthwaite. El Índice de Aridez (IA) se calculó y mapeo para un escenario base (19702000). Posteriormente, se proyectó y mapeó el comportamiento de aridez para algunos escenarios, utilizando los modelos de clima global HADGEM 2.0, GFDLCM 3.0, MIP_ESM y CRNMCM5. Se pronosticaron algunos escenarios, en el mejor escenario la aridez debilitará los ecosistemas húmedos y en el peor escenario aparecerán climas hiper áridos en el área de estudio.
ABSTRACT
An annual mean temperature map was calculated using the Kriging interpolation method for the northcentral zone of Mexico to obtain the current aridity, as well as possible scenarios for the near and distant future. The altitudinal gradient was estimated by linear regression, and it was used to estimate the mean temperature. Climate Influence Areas (CIA) were obtained by superimposing the official precipitation layer and the annual mean temperature layer using Geographic Information Systems tools. Monthly databases of climatic variables were generated for each CIA and potential evapotranspiration was estimated using the Thorthwaite methodology. The Aridity Index (AI) was calculated and mapped for a base scenario (19702000). Subsequently, the aridity behavior of some scenarios was projected and mapped using the global climate models HADGEM 2.0, GFDLCM 3.0, MIP_ESM, and CRNMCM5. Under the best scenario projected, aridity will weaken the humid ecosystems and in the worst scenario, hyperarid climates will appear in the study region.
Keywords: Kriging, Climate Change, Potential Evapotranspiration, Linear Regression.
1. Introduction
The United Nations Framework Convention on Climate Change (UNFCCC) defines climate change as “a change attributed directly or indirectly to human activity that alters the composition of the global atmosphere and affects the natural variability of the observed climate over comparable periods”. In Mexico, aridity is promoted by climate change, which directly impacts population development since it affects the production of goods and services. The advance in aridity can contribute enormously to the loss of land productivity and become a severe problem of desertification. According to data from the United Nations Department of Economic and Social Affairs (ONU, 2018), desertification affects approximately onesixth of the world’s population and 70% of all drylands. Furthermore, desertification decreases food production, increases infertility and soil salinization, increases flooding in lower areas, lead to water scarcity, increases health problems due to windborne dust (e.g., eye infections, respiratory diseases, and allergies), changes biological cycles, and reduces livelihoods, which can contribute to stimulate migration (Hori et al. 2011).
Arid, semiarid, and dry subhumid lands are found in Mexico mainly in the central regions influenced by the Western and Eastern Sierra Madre mountain ranges. Drylands occupy approximately 101.5 million hectares in Mexico, just over half of the territory; arid zones represent 15.7%, semiarid regions 58%, and dry subhumid regions 26.3% (UACh, 2011), so, it is important to study their behavior considering that climate is being modified by humans.
In the present work, aridity is addressed through variables such as temperature, precipitation, and evapotranspiration; statistical tools, geostatistical concepts, model theoretical foundations, theory of spatial autocorrelation, statistical assumptions, and interpolation are described. Different projected scenarios are compared and estimates of current and future weather conditions are presented. The world climate is a whole, and what happens in a particular area may be replicated in other places with similar characteristics; the present study contributes to validate results from other sites worldwide.
2. Study region
Eight economic regions of Mexico are described in reports by the National Commission for the Exploitation and Use of Biodiversity (Bassols, 1985; CONABIO, 2010). Each region can include more than one state that share characteristics of the biophysical environment and with common economic activities in the primary, secondary, and tertiary sectors. In the present work, region II is considered given its current characteristics as dry climates and limiting surface of water distribution, conditions that make it especially vulnerable to climate change from the point of view of drought (CONAGUA, 2021). Region II (Figure 1) covers the states of Chihuahua, Coahuila, Durango, San Luis Potosí and Zacatecas and it is located between 110 º 0 ‘0’ and 95 º 0 ‘0’ west longitude and 35 º 0 ‘0’ and 30 º 0 ‘0’ north latitude.
Fig. 1. Economic regions of Mexico (Bassols, 1985 and CONABIO, 2010).
3. Materials and methods
Precipitation and temperature data were obtained from stations grouped in the Rapid Climate Information Extractor (ERIC) of the National Institute of Water Technology (IMTA) version III (IMTA, 2010). Data quality of many climatological stations is not good, with large gaps of missing data due to little maintenance and even permanent deactivation. Only 496 stations with 80% or more viable information were considered. The analyzed variables were maximum temperature (Tmax) and minimum temperature (Tmin) for the period 19702000. Precipitation values were obtained from the Informatics Unit for Atmospheric and Environmental Sciences (UNIATMOS) of the National Autonomous University of Mexico from 1950 to 2000 (UNAM, 2017).
Additional information was obtained from cartography of climate change models, scenarios of climatology (Hijmans, et al., 2005), climate global models, and emission scenarios of the Intergovernmental Panel on Climate Change (IPCC) that were processed and interpolated for Mexico and Central America, based on the Model for the Assessment of Greenhousegas Inducted Climate Change – MAGICC (UCAR, 2005).
The concept of spatial autocorrelation is supported by Tobler’s principle, which states that in a geographic space everything is related and the closest areas are more related (Tobler, 1970). Spatial autocorrelation helps to understand the degree of similitude of one object with another (ESRI, 2019). In general, there is spatial autocorrelation (SA) whenever there is a pattern in the behavior of the variable according to the geographical location of the data.
Geostatistical techniques require secondorder stationarity data, that is, at least the variance must be the same in different subregions. The lack of stationarity is commonly linked to anomalies in space, or to the existence of a trend or spatial gradient whose dimension is greater than the study region (Cressie, 1993). Stationarity can be a problem when interpolating points in space, but it does not justify abandoning geostatistics in favor of other interpolation techniques such as inverse distance weighted (IDW) since it is equally sensitive to lack of stationarity as described by Isaaks and Srivastava (1989).
The global Moran Index (Moran, 1950) measures spatial autocorrelation based on locations and values of the entities simultaneously. Given a set of entities and an associated attribute, it evaluates whether the expressed pattern is clustered, sparse, or random (ESRI, 2019). The Moran index (I) for spatial autocorrelation is given by equation (1):
(1)
where, z_{i} = (x_{i} – x̅) is the deviation of an attribute or characteristic from its mean; x_{i} is an attribute value for characteristic i; w_{i,j} is a spatial weight between feature i and j, n is the total number of features, and S_{0} is the aggregate of all spatial weights which is computed as follows, . To test if the spatial configuration occurs at random (H_{0}) vs the alternate hypothesis (H_{1}) that the spatial configuration does not occur at random, we use the standardized Moran index, .
The GetisOrd general statistic (G) (Getis and Ord, 1992) measures the concentration of high or low values in a given study region and is given by equation (2):
(2)
where x_{i} and x_{j} are attribute values for characteristics i and j, w_{i,j} is the spatial weight between characteristic i and j, n is the number of characteristics in the data set. The null hypothesis for the High / Low Clustering statistic (General G) states, H_{0 }: There is no spatial grouping of entity values vs the alternate hypothesis H_{1}: There is a spatial grouping of entity values. The test statistic is based on the standardized GetisOrd index given by with . The z_{G} statistic is important for the interpretation when the null hypothesis is rejected. If z_{G} is positive, high values are clustered, otherwise low values are clustered. The software ArcMap was used to calculate Moran’s Index (I) and GetisOrd index (G), the zscores, and the pvalues to assess the significance (ESRI, 2019).
Kriging is an advanced geostatistical procedure that generates an estimated surface from a set of scattered points with zvalues. The Kriging tool performs an interactive investigation of the spatial behavior of the phenomenon represented by the zvalues before selecting the best estimation method to generate an output surface. The kriging technique assumes that data are correlated in space. The kriging method is considered to be BLUE (Best Linear Unbiased Estimator) and is similar to Inverse Distance Weighted. The general equation is given by:
(3)
where: Z(S_{i}) is the value measured at location i; λ_{i} is an unknown weight for the value measured at location i; s_{0} is a location and n is the number of measured values.
In ordinary kriging, λ_{i}, depends on a fitted model, the distance to the location, and the spatial relationships between the measured values around the location. The kriging interpolation method describes dependency rules with variograms and covariance functions (spatial autocorrelation). The kriging methodology assumes that data follow a Gaussian distribution, present stationarity, are spatially continuous, and the autocorrelation is known, only then, kriging is considered an optimal predictor (Gallardo, 2006). In this work, ordinary kriging was chosen due to the nature of the data, in addition, it offers the possibility of making different adjustments through preestablished models.
For realizations of a variable Z at points x_{i}, i = 1,2,… n, say, Z (x_{1}), Z (x_{2}),… Z (x_{n}), we want to estimate Z (x_{0}), and there was no measurement at x_{0}. In this circumstance, the ordinary kriging method estimates the variable as a linear combination of the n random variables as follows:
(4)
where the λ_{i} represent weights of the original values. Ẑ (x_{0}) is the best linear estimator because the weights are obtained in such a way that they minimize the variance of the estimation error; other interpolation methods such as inverse distances or polygonal do not guarantee minimum estimation variance (Samper and Carrera, 1990). Weights are estimated by minimizing V [Ẑ (x_{0}) – Z(x_{0})] subject to , using the method of Lagrange multipliers. The weights λ, also, can be estimated through the semivariance function, for which it is required to know the relationship between the covariogram and semivariance functions. The modeling of semivariograms includes two fundamental stages (Xie and Myers, 1995): estimation of the empirical semivariogram, and fitting of a theoretical model, to determine the descriptive parameters of the semivariogram that will later be used in the estimation (Journel and Huijbregts, 1976).
Temperature varies as a function of altitude; the empirical method developed by Fries et al. (2012) estimates temperature in mountainous regions, calculates the ranges of variation of the annual temperature at different elevations by a simple linear regression model:
(5)
where y_{i} is the temperature (in ºC), x_{i} is the elevation (in meters above the sea level), i = 1,…, n the number of data points, α is the intercept or temperature at sea level, while β is the temperature variation for each meter that increases the elevation, as in nature there is an inverse relationship between elevation and temperature, the value of β is negative, e_{i} is the random error, we assume e_{i}~ NIID (0, σ_{e}^{2}), where “NIID” stands for Normal Independent and Identically distributed and σ_{e}^{2} is the variance associated to the errors. The elevation was estimated with models of geographical elevation provided by the National Institute of Geography and Statistics (INEGI, 2020).
The mean annual temperature layer was obtained with ordinary kriging interpolation; ranges of 1ºC were generated. The Climate Influence Areas (Gómez et al., 2008) were obtained by cross correlating the cartography of mean annual temperature with the estimated annual precipitation. Potential Evapotranspiration (E) was calculated with Thornthwaite method as follows. First, the monthly heat index was calculated from the monthly mean temperature (T_{i}, º C), as ; subsequently, the annual heat index (I) was calculated with equation 6:
(6)
The monthly potential evapotranspiration (mm/month) for months of 30 days and 12 hours of sun (theoretical), E “uncorrected” was calculated also based on temperature data, as:
(7)
where a = 675(10^{–9}) (I^{3}) – 771(10^{–7}) (I^{3}) + 1792 (10^{–5})(I) + 0.49239 is a theoretical adjustment coefficient. Finally, the calculation of the corrected E for month i is obtained through equation 8,
(8)
where, N_{i} is the maximum number of hours of sunshine in month i and d_{i} is the number of days in month i. The Aridity Index (AI) is a qualitative characteristic of the climate, which makes it possible to measure the degree of sufficient or insufficient precipitation. The Aridity Index (AI) was established in in Mexico through the Official Gazette of the Federation, on 1 June 1995 to define an aridity classification based on the Precipitation (P) and Evapotranspiration (E) as AI=P/E (DOF, 1995). The AI allows the classification into eight categories: hyperarid (< 0.051), arid (0.0510.201), semiarid (0.2010.501), dry subhumid (0.5010.651), humid subhumid (0.6510.751), humid (0.7511.251), very humid (1.2512.5) and per humid (> 2.5).
Climate change scenarios are representations of future climate, based on an internally coherent set of climatologic relationships, which are constructed to investigate the potential consequences of climate change. Scenarios were constructed by superimposing layers of annual mean temperature and precipitation; and then databases were generated for the potential evapotranspiration and the AI. The centroids were obtained with automatic tracing of isotherms and isohyets; finally, global climate models were used to project data for near future (20152039) and far future (20752099), considering two values (4.5 and 8.5) of the Representative Concentration Pathway (RCP).
4. Results and discussion
The zscore of the Moran Global Index and G Index of Getis Ord, constructed with a weighting matrix under the Knearest neighbor’s method, were 3.7797 and 5.5738 respectively; in both cases there is less than 1% probability that the pattern is the result of a coincidence, so there is indeed a spatial autocorrelation grouped in clusters for the case of the calculated mean annual temperature.
Fries, et al. (2012) indicate that it is necessary to consider the altitudinal gradient to obtain a correct temperature layer once the interpolation has been carried out, so the fitted linear regression model was:
(9)
The residuals were distributed randomly and with a presumably normal histogram, so we can conclude that the model complies with the assumptions of normality and homoscedasticity. Once the linear regression model was fitted, the equation proposed by Fries, et al. (2012) was used to predict the temperature:
(10)
where, T_{month} corresponds to the mean annual temperature in º C; Γ corresponds to the altitudinal gradient calculated in the regression model (β̂ = Γ = –0.0042); T_{Det} corresponds to a reference altitude (in this case 2000 meters above the sea level, used in the kriging methodology); and Z_{est} corresponds to the elevation of each weather station in meters above the sea level.
The Kriging procedure uses transformations in case data do not comply with the normality assumption, but it is difficult to use transformed values since it will later be necessary to work with the precipitation variable. A ShapiroWilk test was performed to verify data normality, with the hypothesis H_{0}: The population is normally distributed vs. H_{1}: The population is not normally distributed. A test value W = 0.97707 was obtained, with a pvalue < 1 × 10^{–5}, so setting the significance of the test at α = 0.05, we conclude that the population does not come from a normal distribution and a type of transformation is required. To correct the problem of normality, a subdivision of the study region was proposed, considering the physiographic provinces proposed by CONABIO, (CONABIO 2001) as shown in Figure 2.
Fig 2. Physiographic sub regions and distribution of meteorological stations.
The physiographic provinces reduce extreme values and have homogeneous distributions. The ShapiroWilk goodness of fit test is shown in Table I, indicating that none of the subregions has normality problems. Once the assumptions were corroborated, the variograms of each physiographic province were obtained, as well as the adjustment of different models depending on the dispersion of Gamma values, as shown in Table II.
Table I. Results of the Shapiro Wilk goodness of fit test.
Sub Region 
Shapiro Wilk’s W value 
Pvalue for normality test 
R1  0.9893  0.3298 
R2  0.9832  0.1325 
R3  0.9836  0.2689 
R4  0.9824  0.1659 
Table II. Summary of adjustment of variograms, models marked in bold have the best fit and were selected.
Physiographic Region 
Model  SSR  Variogram and Fitted Model  Variogram attributes 
R1 Center 
Exponential Gaussian Spherical Circular 
63.21 60.62 60.16 60.44 
Range:100 Psill: 0.9 Nugget:0.7 

R2 Western Sierra Madre 
Exponential Gaussian Spherical 
38.97 326.57 36.64 
Range:300 Psill: 4.0 Nugget:2.1 

R3 Eastern Sierra Madre 
Exponential Gaussian Spherical Circular 
63.05 88.06 99.30 103.73 
Range:200 Psill: 2.0 Nugget:1.0 

R4 North 
Exponential Gaussian Spherical Circular 
95.43 18.39 18.61 18.11 
Range:100 Psill: 3.2 Nugget:2.1 
Model validation through different error indicators was carried out by two methodologies: Crossvalidation and validation. Table III presents the summary of crossvalidation, obtained through a random sample made up of 80% of the data, the variable temperature was estimated and validated with the stations not included in the sample. In general, values for mean, root mean square, standardized mean, root mean square standardized, and average standard error are small and some of them close to zero, therefore it can be concluded that the adjusted model seems to be correct.
Table III. Summary of the crossvalidation indicators, where R1, R2, R3 and R4 represent the previously defined sub regions.
Topic  R1  R2  R3  R4 
Sample  117  142  173  98 
Mean  –0.0008  –0.0157  0.0366  0.0245 
Root Mean Square (RMS)  0.905  1.5362  0.0258  1.2813 
Standardized mean  –0.0010  –0.0082  0.0209  0.0174 
RMS Standardized  1.0858  1.0823  0.9856  0.8987 
Average Standard Error  0.8319  1.3954  1.1772  1.4239 
Regression function  –0.006*x + 17.17  –0.52 * x + 8.70  0.30 * x + 11.51  0.54 * x + 7.43 
For the validation methodology, the data were divided into 2 subsets: a training set and a test set, the subsets are large enough to be representative and to generate statistically significant results. Again the estimations for average error, root mean square, and average standard error are close to zero, which indicates a good fit (Table IV). Once the validation of the models is completed, the temperatures are interpolated and the altitude gradient is adjusted.
Table IV. Summary of the values of the errors obtained in the validation.
Topic  R1  R2  R3  R4 
Average error  0.075  –0.333  –0.173  0.614 
Root mean square  0.829  1.420  1.193  1.420 
Average standard error  0.089  –0.240  –0.194  0.429 
Subregion R1, mainly presents temperatures between 14 ºC and 19 ºC, with some surfaces between 8 ºC and 10 ºC in the elevated regions, as well as temperatures of 22 ºC to 24 ºC in the northwest. For subregion R2, there are low average temperatures throughout the summit of the western Sierra Madre with temperatures ranging between 8 ºC and 10 ºC, with the highest temperatures also occurring towards the east of the mountain, ranging between 24 and 26 º C. In subregion R3, the highest temperatures are obtained in the southeast, coinciding with the warmest region of the state of San Luis Potosí, with average temperatures of 25 ºC to 26 ºC, as well as the coldest in the center with isotherms of 9 ºC to 10 º C. In subregion R4, temperatures are generally observed ranging between 16 ºC and 18 ºC with the warmest surfaces to the east with temperatures ranging between 20 ºC and 23 º C (Fig. 3).
Fig. 3. Temperature map of region II.
The total annual precipitation was obtained using map algebra
(11)
where i corresponds to the imonth evaluated. The result, reclassified to every 100 mm of rainfall, is presented in Figure 4. Twentyone classes of precipitation were obtained, the most intense colors of blue being the regions with the highest rainfall.
Fig. 4. Precipitation map of region II.
The precipitation is distributed as expected, the most humid regions are those corresponding to the windward part of the western Sierra Madre and the southern part of the state of San Luis Potosí. On the other hand, the regions with less precipitation correspond to the leeward part of both the western and eastern Sierra Madre, coinciding with the theory of Gómez et al. (2008).
Monthly databases of climatic variables were generated for each CIA. Centroids of each polygon were obtained using geometric methods. The aridity index was calculated with estimations of annual mean temperature (Tmed), potential evapotranspiration (E), and precipitation (P). The baseline scenario was built with climate data from the 19702000 period. The data for the period 20002020 were not used because it shows a considerable decrease in quality in most of the stations considered. The calculated aridity baseline scenario is presented in Figure 5. The percentages for each classification were as follows: hyperarid (0%), arid (1%), semiarid (54%), dry subhumid (15%), humid subhumid (9%), humid (15%), very humid (6%) and per humid (0%). The baseline scenario for region 2 presents a predominance of a semiarid climate with more than half of the surface region (in km^{2}), the percentages of Perhumid and Hyperarid climates are close to 0%. The highest concentration of humid regions is concentrated in the western highlands, more specifically in the windward region, and the driest regions are in the northcentral region.
Fig. 5. Aridity Base Scenario of region II (19702000).
Table V and Figure 6 show the estimates of four global climate models considering two Representative Concentration Pathway (RCP) and two horizons (H1: 20152039 and H2: 20752099), in general, there is a transition from left to right, from humid classes to arid classes, the most optimistic scenario is CRNMCM5 where the “hyperarid” class is not reached and in the worst scenario, the closest horizon only occupies 6245.40 km^{2} of the “arid” class. The most pessimistic models are MPIESP and GFDLCH3, since in both cases the “hyperarid” class is reached with 5639.26 km^{2} (0.9%) and 5621.52 km^{2} (0.9%) respectively. The class “perhumid” represented by just 1 km^{2} only exists in the base scenario and is extinguished in the rest of the scenarios. A model that is neither very pessimistic nor very extreme is HADGEM 2.0, which predicts an “arid” surface of 349325.99 km2 (53.8%) without reaching the “hyperarid” class.
Table V. Results of estimates in percentage for different horizons in eight climate classes. Four models are presented.
Classification  
Model  RCP  H  Per humid 
Very wed 
Wed  Sub humid Wed 
Sub humid Dry 
Semi arid 
Arid  Hyper arid 
base scenario    0  0.0%  6.3%  14.6%  9.5%  15.1%  53.5%  1.0%  0.0% 
HADGEM 2.0  4.5  1  0.0%  4.8%  14.1%  2.0%  13.5%  64.7%  1.0%  0.0% 
2  0.0%  2.8%  15.6%  2.3%  5.1%  73.1%  1.0%  0.0%  
8.5  1  0.0%  7.6%  12.9%  3.4%  13.7%  61.5%  1.0%  0.0%  
2  0.0%  1.9%  7.0%  1.0%  2.1%  34.3%  53.8%  0.0%  
CRNMCM5  4.5  1  0.0%  9.9%  12.6%  6.4%  9.1%  61.1%  0.9%  0.0% 
2  0.0%  3.6%  16.4%  1.1%  11.3%  66.7%  1.0%  0.0%  
8.5  1  0.0%  9.6%  11.4%  4.6%  12.6%  60.9%  0.9%  0.0%  
2  0.0%  2.8%  7.3%  10.5%  8.8%  69.6%  1.0%  0.0%  
GFDLCH3  4.5  1  0.0%  8.8%  15.4%  5.6%  8.3%  61.1%  0.9%  0.0% 
2  0.0%  8.6%  15.3%  5.0%  9.1%  61.2%  0.9%  0.0%  
8.5  1  0.0%  3.0%  15.4%  1.1%  5.8%  72.9%  1.7%  0.0%  
2  0.0%  2.8%  5.6%  8.8%  6.6%  20.9%  54.4%  0.9%  
MPIESP  4.5  1  0.0%  7.3%  12.9%  2.6%  9.7%  66.5%  1.0%  0.0% 
2  0.0%  3.5%  7.0%  10.4%  9.9%  68.3%  1.0%  0.0%  
8.5  1  0.0%  3.6%  16.8%  2.4%  9.1%  67.2%  0.9%  0.0%  
2  0.0%  1.9%  7.9%  0.3%  10.5%  24.4%  54.2%  0.9% 
H: Horizons, 0: 19702000; 1: 20152039; 2: 20752099.
RCP: Representative Concentration Pathway.
HADGEM2.0: Met Office Hadley Centre (MOHC) (United Kingdom model).
CRHMCM5: Centre National de Recherches Météorologiques (French model).
GFDLCH3: Geophysical Fluid Dynamics Laboratory (USA model).
MPIESP: Max Planck Institute for Meteorology (Germany model).
Fig. 6. Mapping of aridity evolution scenarios for models: a) HADGEM2.0; b) MIPESP; c) CRNMCM5; d) GFDLCM3; I: RCP=4.5, Horizon 20152039; II: RCP=4.5, Horizon 20752099; III: RCP=8.5, Horizon 20152039; IV: RCP=8.5, Horizon 20752099.
In summary, the results indicate that all the assumptions of the Kriging model were fulfilled and, therefore, the best unbiased linear estimator was developed. The validation model showed that estimates were correct. A base aridity index and projections from global climate models were obtained for near and far futures. Such projections are useful to make decisions in the present and change the projected trend. All the global climate models predict an increase in aridity for the far future, and RCP = 8.5, so we should take action to avoid large RCP levels. The information generated in this study has technical foundations.
For future research, the study can be extended to develop projections for all regions of Mexico. Having reliable information helps us make the best decisions and, in terms of climate, it allows us to reduce climate change effects.
5. Final remarks
The availability of climate information in Mexico is severely affected by different factors, such as meteorological stations with no maintenance and less than 70% of the data. The aridity in regions where land has severe degradation problems can lead to a loss of productive capacity practically irreversible. The use of geostatistical interpolation methodologies for unsampled points is very useful to model the behavior of different climatic variables. Temperature and precipitation detailed maps were obtained for economic region II of Mexico. Different scenarios were predicted with global climate models. Results show adverse scenarios: in the best prediction, aridity will strongly weaken the humid ecosystems, and in the worst scenario hyperarid climates will appear, practically inhospitable for many of the species with limited environmental resilience that today inhabit these regions. It is urgent to carry out specialized studies of territorial order, natural resources management, regeneration, and adaptation work for the regions where ecosystems are most vulnerable. Some proposals that can be highlighted are: wellmanaged reforestation plans, constant cooperation with multilateral organizations, transition towards clean and renewable energies, reducing emissions into the atmosphere, and reducing the potential of the highest Representative Concentration Pathways.
Although Mexico does not have an atmospheric circulation model, some efforts have been made to represent the country’s environmental conditions, and this study, supported by geostatistics, is one of them.
References
Bassols BA. 1985. Recursos naturales de México: teoría, conocimiento y uso. Editorial Nuestro Tiempo. Segunda Ed.
CONABIO. 2001. Provincias Fisiográficas de México, México DF.
CONABIO. 2010. Catálogo de metadatos geográficos. Comisión Nacional para el Conocimiento y Uso de la Biodiversidad. Available at: http://www.conabio.gob.mx/informacion/metadata/gis/recomgw.xml?_xsl=/db/metadata/xsl/fgdc_html.xsl&_indent=no (accessed 2021 January 3)
CONAGUA. 2021. Monitor de Sequía en México. Available at: https://smn.conagua.gob.mx/es/climatologia/monitordesequia/monitordesequiaenmexico (access 2021 May 19)
Cressie N. 1993. Statistics for Spatial Data. Revised Edition. John Wiley and Sons. New York.
DOF. 1995. Diario Oficial de la Federación, (DOF: 07/06/1995). Available at: http://dof.gob.mx/nota_detalle.php?codigo=4875100&fecha=07/06/1995 (accessed 2021 May 25)
ESRI. 2019. Cómo funciona Autocorrelación espacial. Available at: https://desktop.arcgis.com/es/arcmap/latest/tools/spatialstatisticstoolbox/hhowspatialautocorrelationmoransispatialst.htm (accessed 2021 January 3)
Fries A, Rollenbeck R, Naub T. 2012. Near surface air humidity in a megadiverse Andean mountain ecosystem of southern Ecuador and its regionalization. Agricultural and Forest Meteorology, 152: 1730. https://doi.org/10.1016/j.agrformet.2011.08.004
Gallardo A. 2006. Geostadística. Ecosistemas, 15(3): 4858. https://www.revistaecosistemas.net/index.php/ecosistemas/article/view/161
Getis A, Ord JK. 1992. The Analysis of Spatial Association by Use of Distance Statistics. Geographical Analysis , 24: 189206. https://doi.org/10.1111/j.15384632.1992.tb00261.x
Gómez JD, MonterrosoRivas A, Gay C. 2008. Spatial estimation of mean temperature and precipitation in regions of scarce meteorological information. Atmósfera, 21: 3556. https://www.revistascca.unam.mx/atm/index.php/atm/article/view/8596
Hijmans RJ, Cameron SE, Parra JL, Jones PG. 2005. Very high resolution interpolated climate surfaces for global land areas. International Journal of climatology, 25: 19651978. https://doi.org/10.1002/joc.1276
Hori Y, Stuhlberger C, Simonett O. 2011. Desertification: a visual synthesis. United Nations Convention to Combat Desertification (UNCCD), First edition, GRAPHI 4 Press, Bresson France.
IMTA. 2010. Extractor Rápido de Información Climatológica, México: Instituto Mexicano de Tecnología del Agua.
INEGI. 2020. Continuo de Elevaciones Mexicano (CEM). Available at: https://www.inegi.org.mx/app/geo2/elevacionesmex/ (accessed 2021 January 3).
Isaaks E, Srivastava R. 1989. An introduction to applied geostatistics. Oxford University Press.
Journel AG, Huijbregts CJ. 1976. Mining Geostatistics. United Kingdom.
Moran PAP. 1950. Notes on Continuous Stochastic Phenomena. Biometrika, 37: 17–23. https://doi.org/10.2307/2332142
ONU. 2018. Ordenación de los ecosistemas frágiles: Lucha contra la desertificación y la sequía. Departamento de Asuntos Económicos y Sociales. Available at: https://www.un.org/spanish/esa/sustdev/agenda21/agenda21spchapter12.htm (accessed 2021 January 3)
Samper CFJ, Carrera RJ. 1990. Geoestadística. Barcelona: Centro Internacional de Métodos Numéricos en Ingeniería.
Tobler WR. 1970. A computer movie simulating urban growth in the Detroit region. Economic Geography, 46: 234240. https://doi.org/10.2307/143141
UACh. 2011. Actualización y delimitación de zonas áridas, semiáridas y subhúmedas de México. Reporte final de proyecto de investigación. Departamento de Suelos. Universidad Autónoma Chapingo.
UCAR. 2005. Model for the Assessment of Greenhousegas Inducted Climate Change. A Regional Climate Scenario Generator. Available at: https://www.cgd.ucar.edu/cas/wigley/magicc/ (accessed 2021 January 3)
UNAM. 2017. Portal de datos de escenarios de cambio climático. Responsorio Institucional Instituto de Ciencias de la Atmósfera y Cambio Climático, UNAM. Available at: https://ri.atmosfera.unam.mx:8150/geonetwork/srv/eng/catalog.search#/home (accessed 2021 January 3)
Xie T, Myers DE. 1995. Fitting MatrixValued Variogram Models by Simultaneous Diagonalization (Part I: Theory). Mathematical Geology, 27(7): 867875. https://doi.org/10.1007/BF02087100