Proportional Variation of Potential Groundwater Recharge as a Result of Climate Change and Land-Use: A Study Case in Mexico

This work proposes a methodology whereby the selection of hydrologic and land-use cover change (LUCC) models allows an assessment of the proportional variation in potential groundwater recharge (PGR) due to both land-use cover change (LUCC) and some climate change scenarios for 2050. The simulation of PGR was made through a distributed model, based on empirical methods and the forecasting of LUCC stemming from a supervised classification with remote sensing techniques, both inside a Geographic Information System. Once the supervised classification was made, a Markov-based model was developed to predict LUCC to 2050. The method was applied in Acapulco, an important tourism center for Mexico. From 1986 to 2017, the urban area increased 5%, and by 2050 was predicted to cover 16%. In this period, a loss of 7 million m3 of PGR was assumed to be caused by the estimated LUCC. From 2017 to 2050, this loss is expected to increase between 73 and 273 million m3 depending on the considered climate change scenario, which is the equivalent amount necessary for satisfying the water needs of 6 million inhabitants. Therefore, modeling the variation in groundwater recharge can be an important tool for identifying water vulnerability, through both climate and land-use change.


Introduction
Land-use cover change (LUCC), especially the transformation of natural vegetation to agricultural and urban uses, can directly modify the hydrological cycle, resulting in increased flooding and decreased groundwater recharge, among other effects [1][2][3][4]. In agricultural areas, the topsoil compaction leads to a loss of the water retention and infiltration capacity of the subsurface. In urban areas, the introduction of relative impermeable materials or surfaces prevents groundwater recharge [5,6].
Precipitation, lateral flow from adjacent unconfined aquifers, irrigation practices and pipeline leakage are the main sources of groundwater recharge [7][8][9]. The vertical water that drains from the surface to the root zone (unsaturated zone) is called infiltration [10], and the subsequent displacement of water by capillarity and gravity to the groundwater (saturated zone) is called percolation [9]. Hence, for this study, the potential groundwater recharge (PGR) is related to the vertical flow from precipitation that has the potential to percolate to the aquifer [8,9], and the main factors that allow topography of the region is rugged, ranging from 0 to 2079 m.a.s.l. [50]. The main water bodies include the Sabana and Papagayo Rivers, Tres Palos Lagoon, and the El Veladero Protected Natural Area, where human settlements are prohibited ( Figure 1) [51]. Accordingly, this area has not experienced LUCC, and was considered as restricted by LUCC in the analysis. This region covers two aquifers, Acapulco Bay and La Sabana. Both present similar characteristics, as follows: unconfined aquifers, hydraulic conductivity from 0.12 to 9.12 × 10 −3 m/s, and transmissivity between 1.5 and 120 × 10 −3 m 2 /s. About lithology, the aquifers show a predominance of unconsolidated sediments and metamorphic rocks, and according to official sources, the information on depth to the water table is scarce. Some data indicate values from just a few centimeters in the bay to 9 m at the urban zone [52,53]. With respect to tourism, the main infrastructure is concentrated around the Bays of Acapulco and Puerto Marqués Bay. In 1950, 246 ha were identified as belonging to the tourism sector by Carrascal and Pérez [54]. Since then, this area had increased to 450 ha for 1989. From 1987 to 1993, the most significant tourism development of Acapulco was built: The Zona Diamante. This urban project involved a land-use change of 270 ha, extending from Puerto Marqués Bay toward the western portion of the municipality between the Pacific Ocean and Tres Palos Lagoon [54]. Several other tourism and residential buildings have been constructed in low-lying areas near the coast. The most susceptible areas to flooding are located in the southern portion of the municipality and to the southwest of Tres Palos Lagoon. Notably, two large hotel complexes with an extension of approximately 133.7 ha along Puerto Marqués Bay are located at no more than 20 m.a.s.l. and, therefore, are highly susceptible to flooding [44,54].
Concerning some explanatory variables of LUCC (Table 1), roadways are triggers of urban and agricultural growth [55,56] because they enable higher access and the transport of merchandise. Existing human settlements tend to expand given their proximity to existing public services, among other factors [57]. For this reason, buffers based on the distance to urban centers in 2003 and roadways were created [48,58]. Given the rugged topography of the area, it is also important to consider the slope and elevation as explanatory variables. Given that roadways reinforce urban LUCC because they provide access to services and businesses, the national roadway networks of 2014-2017 were also considered as an explanatory variable for LUCC.

Materials and Methods
This section ( Figure 2) presents the supervised classification methodology for determining landuses from 1986, 2003 and 2017, as well as the determination of LUCC in 2017 using a predictive model, and its validation of current land-uses. The model was then used to forecast land-uses by 2050. Finally, the variation in potential groundwater recharge during the study period was estimated through a distributed model. .  Table 2 shows the characteristics of the satellite images of Landsat 5 TM, 7 ETM + and 8 OLI (Path/Row 26/48) used in this study for April 1986, 2003 and 2017 respectively, as well as those downloaded from the website of the United States Geological Survey [62], which provides free-access images with adequate geometric correction for multi-temporal analyses [2,61,63]. The year 1986 was selected as it is the year before the construction of the largest urban increases registered in the study area (Zona Diamante), mainly driven by the high tourist demand. Besides, a period of 31 years is suggested as the minimum time required to validate a predictive model of LUCC [64,65]. The spectral bands were segmented for those three years using a watershed-based algorithm in the TerrSet (Version 18.31) software [60,66], and the training sites were chosen in the segmented images. For the first two years, these were selected from well-known sites identified on orthophotos (Table 3). For 2017, they were selected through a field review and, by consulting Google Earth [67,68] following a similar spectral response pattern, visually homogenous and evenly distributed across the images according to Mondal et al. [68]. The training sites were evaluated by estimating the separability of the spectral signatures via the transformed divergence, which is a statistical method that compares the mean and co-variance of pairs of spectral signatures to evaluate the accumulated difference between the training sites of the selected land-uses. Accordingly, it was possible to perform the supervised classification using the maximum likelihood algorithm [71,72].

Land-Use Change Forecasting
The land-use change forecasting consisted of the selection and evaluation of the explanatory variables of LUCC, such as the urban infrastructure, topography and climate. Liang et al. [56] mentioned that the distance to roadways and human settlements is a crucial determinant of LUCC [61]; also, orographic characteristics, such as slope and elevation, define the speed of agricultural and urban expansion [73]. To further evaluate the association between the LUCC and explanatory variables, the Cramer's V test was used, considering values between 0.15 and 0.4 as indicators of an acceptable relationship between the land-use change and the explanatory variable [58,66,73].
A neuronal network algorithm was applied to calculate the transition potentials. This algorithm is one of the most applied because it can solve complex variables and provide highly precise results [74,75]. The forecasting process was performed using Markov chains [61, 74,76]. The land-uses in 2017 were forecasted and validated using the Kappa statistic [74,77] and the ROC curve, which compared the coincidence between the forecasted and actual land-uses [77]. Once we had validated this model, an image of land-use in 2050 was constructed.

Estimation of PGR
Estimating PGR through a hydrological balance leads to a series of decisions, which must be made on the following: the selection of methods for determining evapotranspiration, runoff, etc.; hydrological models providing information with adequate space and time formats; and even the complexity of the model. As this study focuses on LUCC, it requires a model able to provide information with a territorial approach that is as simple as possible, to keep the PGR projection in the future lying just within the forecast of temperature, rainfall and LUCC.
According to its parameters as a function of space, a model could be semi-distributed or distributed. For instance, it is possible to find semi-distributed models such as HBV-model [78], the Variable Infiltration Capacity (VIC) model, which is recommended for humid and agricultural areas [79,80], or the soil and water assessment tool (SWAT), which focuses on agricultural areas [80]. There is also TOPMODEL, which is used for moderate topography [80,81].
These models are often calibrated and validated in a given time, and when they are applied throughout long periods, the land-use is assumed to be time-invariant. In this case, it is implausible to assume that the model parameters are still suitable when large extensions of land-use have changed [23]. Therefore, in the present study a distributed model is preferred, where parameters, inputs and outputs can vary spatially.
Among the distributed models and according to Devi et al. [80], the conceptual ones are simple and can be easily implemented in computer code, unlike the psychically-based ones that are complex, require human expertise and computational capability, and suffer from scale-related problems.
The conceptual distributed models are usually composed of a production function (or generation routine) and a transfer function (or concentration routine). This approach allows one in the first place to estimate vertical flows (e.g., PGR) in each cell of the distributed model grid, independent from characteristics changes (such as LUC) on the other cells. Secondly, the transfer function estimates the runoff concentration as an output that could be validated with field measures.
Some models work inside a Geographic Information System (GIS) environment; for example, InVEST Water Yield Model is a model that estimates the water yield as a percentage of precipitation based on the Budyko curve [82]. Nevertheless, no distinction is made between a sub-surface and surface flow [83]. Given this lack of distinction, in this work, the production function for PGR estimation (step 6) and its variation over time (step 7) is carried out by another GIS environment module [84]. As well as other models, this tool estimates PGR by abstracting from the precipitation, the potential evapotranspiration (PET) and the surface runoff monthly. Precipitation imagery every month from 1980 to 2018 was obtained from the databases of Thornton et al. [85].
As regards PET, the Thornthwaite method [86] was considered adequate due to its impact on runoff assessment under climate change and its feasibility in the study zone. Variation in rainfall and temperature as a result of climate change modifies evapotranspiration [87] and the overall hydrological cycle. Under climate change, Yates and Strzepeck [88] found this method to be betterranked compared to Hargreaves', and just below the Priestley-Taylor method in terms of the standard error. The latter one, as well as the Penman method, is a physically-based model, and therefore dataintensive and data-sensitive. Furthermore, the temperature-based methods may be sufficient for large scales, if it is assumed that other climate variables covary with temperature [89]. Pan's et al. [90] calculations demonstrate the global increase tendency of evapotranspiration on the order of 60% to 80% in some regions up to 2090. It is also noticed that different climatological regions show different trends in PET estimations depending on the method [88]. For instance, Thornthwaite tends to underestimate PET in sub-tropics, and this is most likely because the method was originally calibrated in humid environments [89]. This allows for inferring that the method could be feasible in regions with enough relative humidity like the present study zone (> 75%). In this case, Van der Schrier et al. [91] show a difference from -0.1 to 0.1 m/year of PET between Thornthwaite and Penman methods.
On the other hand, the surface runoff is depicted by a coefficient Keij as the ratio between the runoff and rain [92,93], which involves slope, altitude, soil type and land-use.
The soil type was associated with its permeability based on the soil classification of INEGI [94]. The Regosol soil type was linked with high permeability, while the Solonch and Phaeozem soil types were associated with low and medium permeability, respectively. The Ke ranges from 0.8 to 0.35 for low, 0.7 to 0.25 for medium and 0.5 to 0.05 for high permeability depending on slope and land-use [95,96]. A digital elevation model (DEM) with 13.5 m of spatial resolution [50] was used to estimate the slope of the terrain. Later, the slope was categorized into high (>50%), moderate (5% to 20%) or low (0% to 1% [93,97,98]). About the land-use, it was classified according to categories usually found in the literature [8,11,82,99] and following the classification for LUCC: (a) urban zone (populated areas); (b) croplands; (c) forest areas; and (d) water bodies.
The surface runoff was validated throughout the efficiency coefficient Nash-Sutcliffe [100] with available information from gauging stations between 1981 and 1998 [101]. In addition to the runoff generation, the groundwater reaction (baseflow) and human activities influence river discharges. Therefore, given the monthly basis, the baseflow was separated based on separating techniques from the time series of stream flows [102]. Due to the inclusion of human activities that can improve the simulation of monthly discharges according to Veldkamp et al. [103], a monthly discharge from the urban use was added. This discharge was a function of the population per year, the rate of population increase (11,000 inhabitants/year [32]), and a discharge rate of 75% of water supply (90 L/inhabitant/day).
To forecast PGR in 2050, the scenarios proposal regards (a) current annual climate conditions, (b) moderate concentrations of greenhouse gases, and (c) extreme concentrations of greenhouse gases. The second and third scenarios are associated with variations in rainfall and temperature according to the RCP 4.5 and RCP 8.5 models of climate change [104]. These scenarios were estimated using the Max Planck Institute Earth System Model (MPI-ESM) [49,105].
In step 7, the groundwater recharge volumes were estimated for each year (1986,2003,2017, and 2050). The variation in potential groundwater recharge ∆PGRT-t (Equation 1) between years T and t was determined as a function of the potential groundwater recharge volume PGRij of the pixel i in month j.
(1) Figure 3a shows a true-color image of the study zone where one can see a wide range of spectral signature as a function of the cover land. However, another type of band composition allowed for detecting specific land covers in a more differentiated manner (phase 1). For example, the false-color composition 234 (Figure 3b) highlights the differences between vegetation and urban area, as well as among several types of vegetation. Figure 3c shows a 345 composite image that was useful to differentiate between agricultural and natural vegetation [60]. All the selected satellite images satisfied a cloud cover lower than 20%. In addition to composite images, other materials were necessary in order to define training sites for supervised classification. The orthophotos obtained for 1995 and 2007 provided certainty about the land cover in past times because of its higher resolution concerning those from the satellite. For example, Figure 3d,e show how the construction of urban infrastructure followed the pattern drawn for different types of vegetation found in the composite images of 1986.

Supervised Classification
Based on the spectral signatures of the training sites, an analysis of separability was carried out using the transformed divergences, resulting in separability values of 1829 to 2000 units. These values are acceptable according to Richards and Jia [106]. The lowest separability value (1829) belongs to the urban/bare soil areas.
Urban growth is mostly observed in the southwestern portion of the municipality bordering the Pacific Ocean and toward the southeastern portion where the El Veladero Protected Natural Area is located. Further, urban growth appears to follow the pathway of the main roadways that connect Acapulco with central Mexico (Figure 4). Overall, most of the LUCC is associated with an increase in the urban area (Table 4). From 1986 to 2017, Acapulco expanded its urban zone to more than double its original size, mostly at the expense of agricultural areas. Further, the forest area decreased by 2%. Some forest areas were displaced by agriculture, although not a significant amount. The overall decrease in agricultural areas might be attributed to the low agricultural productivity of the region in comparison to the tourism sector. These changes imply an increase in the impermeable surface. However, no perceptible changes to the water bodies of the Papagayo River and Tres Palos Lagoon can be distinguished.

Forecasting of Land-Use Change and Validation
After the land-uses were determined for 1986, 2003, and 2017, the relationships between the LUCC and potentially explanatory variables were evaluated according to Cramer's V test. The urban infrastructure, climate and orographic variables were strongly related to LUCC, with V values greater than 0.28 and mostly greater than 0.4 ( Table 5). The slope had the strongest relationship with landuse change, followed by the temperature. To generate a predictive model of LUCC, which was first verified for 2017, a probabilistic Markov chain approach was used considering the land-use classifications in 1986 and 2003. In Figure  4, similarity can be observed between the current land-use map of 2017 and the forecasted land-uses by the model. The growth rate from 1986 to 2017 was about 2.88 km 2 /year for the urban zone, whilst the decrease rates for forest and crops were 0.78 and 2.16 km 2 /year, respectively ( Table 6). The underestimated 3% in LUCC for agricultural use (the northern portion in Figure 4) is likely because the crop area was relatively constant from 2003 to 2017. In the validation of the predictive model of land-uses for 2017, an area under the ROC curve (AUC) of 86% was observed, which was higher than the 71% minimum recommended by Pontius [77]. Furthermore, the Kappa value of 0.85 confirmed that the model adequately predicted the change in the areas of land-use [77].
Given the satisfactory results of the model, LUCC was also estimated for 2050. Table 6 shows that the growth rate of urban use between 2017 and 2050 (2.5 km 2 /year) remains slightly below that of the former period. The existing urban area gets expanded toward the northern and northwestern part of Acapulco, mainly at the expense of agricultural areas and, to a lesser extent, forest areas.
Urban development seems to follow a north-south trend along the federal highway to central Mexico. New urban settlements also appear near water bodies, mainly between the Papagayo River and Tres Palos Lagoon. Along the coastline, a thin black line of urban settlements can be observed to extend from one extreme to the other.
The tendencies in LUCC for 2050 are similar to those of the previously analyzed 31-year period. From 1986 to 2017, the forest area decreased by 2% and the agricultural area by 4%. Meanwhile, the urban area grew by 5%. According to the predictions, from 2017 to 2050, the forest area will lose more than 1% (around 24 km 2 ), and the agricultural area will lose 5%. Notably, the urban area will gain 7%.

Current and Future Potential Groundwater Recharge
As is mentioned in the Method, the conceptual and distributed hydrologic model must be validated with field measures. In this context, the validation of surface runoff through the information provided by gauge stations would allow for assuming the coherence among the rest of the hydrologic flows. Figure 5a shows a comparative view of the hydrographs with the simulated and observed surface monthly runoff volume corresponding to the contributing catchment area of gauge station 19,018 (see Figure 1). According to several criteria [105], the goodness of fit for the model is satisfactory because it provides an adjusted coefficient of determination equal to 0.768, and a Nash-Sutcliffe coefficient of 0.767. In annual terms, the average simulated runoff volume was equal to 172.7 hm 3 /year (with a standard deviation of 57.2 hm 3 /year) whilst the average observed runoff volume was equal to 173.5 hm 3 /year (with a standard deviation of 90.0 hm 3 /year). The estimated annual runoff volume for 1986, 2003 and 2017 was 161.61 hm 3 , 172.1 hm 3 and 175.28 hm 3 , respectively. Although the accuracy of the estimated runoff could be improved by increasing the number of elements classified relative to land cover (for instance, different types of forest vegetation instead of just one), it is relevant to regard that the more classes that are included in the supervised classification, the smaller the difference between the spectral signatures of classes could be. , which depicts between 10% and 36% of the annual precipitation volume. The official information reports a proportional recharge of 18% approximately [52,53], and remains in the estimated range. In comparison with other studies, Moeck et al. [107] provide the relationship between the recharge rate and several variables around the globe. One of these variables, the Aridity Index (ETp/P), shows that recharge rates > 500 mm/year from more than 5000 locations are associated with Aridity Index values between 0.8 and 1.4, approximately. In the present study case, the mean recharge rate of 740 mm/year is coherent with the mean Aridity Index resulting from the model, equal to 1.26.
Visually, few significant changes in the spatial distribution of PGR are observed ( Figure 6) from 1986 to 2003. Nevertheless, in 2017 the PGR in the northern and northeastern parts of the municipality shows a clear decrease as a result of the lost forest area (41.6 km 2 ). The most significant decrease in potential groundwater recharge over the initial 31-year period is located around the Tres Palos Lagoon and near the Pacific Ocean. This may be attributed to the transition of several highly permeable areas (forest and agricultural uses) to urban uses. The estimated PGR for 2050, given current climate conditions (Figure 7a), decreased mainly due to the expansion of the urban settlements bordering the Papagayo River and toward the extreme northern part of the municipality. Under these conditions, from 2017 to 2050, PGR is expected to decrease considerably by 65.97 hm 3 , or almost 10 times compared to the 1986-2017 period (Table 7). Under the RCP 4.5 and 8.5 climate change scenarios, the PGR estimated for 2050 decreases even more dramatically across more than half of the municipality (Figure 7b,c). This can mainly be attributed to the forecasted changes in temperature and rainfall (Table 8). Although the two climate scenarios show similar behavior, the RCP 8.5 scenario shows a greater loss of groundwater recharge, corresponding with an annual decrease of more than 100 hm 3 , mostly as a result of the annual decrease in rainfall depth of 150 mm in the southern part of the municipality. With respect to the temperature, a significant increase of about 3° C is observed across nearly half of the municipality in both scenarios, compared to 1970-2000. This increase would result in greater evapotranspiration and, consequently, decreased groundwater recharge across more than half of the municipality.  The loss of groundwater recharge from 1986 to 2003 is only one-third of the loss from 2003 to 2017. The loss from 2003 to 2017 is mostly due to the 3% decrease in the forest area. However, the decrease in groundwater recharge from 2017 to 2050 (241 hm 3 ) is estimated to be 9 to 37 times greater than that from 1986 to 2017. This can be largely attributed to changing climate conditions and, secondly, to urban expansion over areas of high permeability. The contribution of LUCC to the reduced potential groundwater recharge is between 25% and 27%.
Over the analyzed period, the largest infiltration area is associated with forest areas in the northeastern part of the municipality. Groundwater recharge generally decreases in the direction of the Pacific Ocean. Further, a decrease in groundwater recharge is observed in the lower part of Acapulco Bay, along the coastline, and around Tres Palos Lagoon and Puerto Marqués Bay.
Recapitulating the effect of LUCC on the PGR, there still are some elements to analyze. Considering that the urban use of land constituted 5% of the total ground surface from 1986 to 2017 and 7% from 2017 to 2050, this seems to be an insufficient change to account for a PGR variation of almost 10 times in the same periods. However, an evaluation regarding both the change of other land uses and the location of those changes is necessary. Regarding the term elasticity, defined as the rate of the proportional variation in time i of the PGR ∆ , related to the proportional variation in time i of LUCC for forest and crops ∆( + ) , ( + ) , it is possible to appreciate a change of 24% in PGR for each unit of variation in LUCC from 1986 to 2017. In contrast, from 2017 to 2050, the elasticity is 170%, which also denotes a sensitivity in LUCC associated with the loss of both crops and forests.
On the other hand, it was estimated that the loss of PGR (assuming the PGR image of 1986) corresponds to the LUCC associated with urban use for 2017 and 2050. The 89.5 km 2 of urban-use land increased in 2017, leading to a loss of 23 hm 3 /km 2 in PGR, whilst the 115.4 km 2 increase in urbanuse land in 2050 (in relation with 2017) leads to a loss of 29 hm 3 /km 2 in PGR. That means the LUCC matters, but also where this change is generated. In this context, Figure 8 shows the urban LUCC between 1986, 2017 and 2050 as related to the PGR relative to 1986, classified in quartiles. From the urban LUCC in 2017, 7 of the 89.56 km 2 (8%) are located in the 4th quartile (Q4) of the PGR relative to the 1986 image (375-685 mm/year). For 2050, the urban LUCC located in the same quartile reaches 15 km 2 (13% of the increase). The information before-mentioned leads to a series of possible strategies in sustainable water resource management. For instance, in the case of the urban LUCC in the quartile of higher PGR (Q4), it is possible to establish a Payment for Hydrological Environmental Services (PHES) as compensation to preserve areas of hydrologic importance. The PHES has been implemented in Ecuador, Nicaragua, and Mexico, e.g., the government of Mexico paid landowners 2109 USD/km 2 /year for hydrologic conservation areas in 2019 [109]. A second strategy refers to using the PGR image classified in quartiles as a weighted input of a Multicriteria Analysis to minimize the urban land-use growth in the context of Urban Planning and Management [110][111][112][113]. The proposal of these strategies becomes more relevant every time as the population, and consequently its water requirements, grow immeasurably (by 2050, the population growth in Acapulco will reach just over 11,000 inhabitants/year).

Conclusions
The present work analyzed not just the effect of the LUCC on the potential groundwater recharge (PGR), but also the combined effect regarding climate change. Although it is suggested to assess this kind of effect as a function of the location, altitude, and other climatic characteristics, as expected for a sub-tropical region, the climate change showed a higher proportional impact than the LUCC.
One of the modeling keys for projecting the impact on PGR throughout the period was keeping the driver variables as simple as possible in the forecasting of both LUCC and PGR estimation. On one hand, the classification of land cover into a few groups allowed the elements in each group to be more similar in their potential capacity for infiltration, whereas the difference from other groups becomes higher. On the other hand, a distributed model for PGR estimation is kept independent of the vertical flows among the cells, so the global value is the sum of every flow instead of the result of a combination of attributes in a whole area. Besides, the use of empirical rather than physically-based models to obtain evapotranspiration meant the PGR could be estimated by variables associated with climate change projection, assuming the rest of the variables covary with temperature. Otherwise, the model would increase the complexity of inferring the variation in other physical variables, such as specific humidity, among others.
For the assessment of the above-mentioned effects, a platform in a Geographic Information System environment provided a powerful tool in the same way, as it requires the variables to be in a distributed format. In this case, remote sensing was demonstrated to be an adequate source for forecasting LUCC over time. Nevertheless, it still needs specific assistance in order to identify accurately training sites. Besides the knowledge and expertise of the researcher concerning the study zone, other complementary materials supported this task: orthophotos for older periods and new developments, and Google Earth for recent ones.
One notable finding is the complex dynamics of LUCC over time. For instance, and opposite to what was expected in the study zone, the forest area slightly increased over the first period (1986 to 2003) but decreased in the second period (2003 to 2017). The agricultural area slightly increased in the second period after decreasing in the first period.
In addition to the classified images, the LUCC required some explanatory variables to generate the forecasting. The one with more influence in the study area was the slope, owing to rugged geography. Urban growth does not occur entirely in areas with higher slopes. Besides, high altitude, low temperature and high precipitation converge in the study zone with low settlements. This analysis could be improved by considering additional socioeconomic variables (for example, the price per m 2 of residential and marginal areas), as suggested by Yalew et al. [16]. However, the satisfactory LUCC forecasting showed that its direction and magnitude have direct and indirect socio-economic implications, especially considering the constant growth of the population and tourism sector in Acapulco de Juárez.
Once the LUCC was successfully validated through a Cramer's test and the area under the curve ROC was determined, the applied distributed model to estimate the potential groundwater recharge provided satisfactory results, according to observed measures.
According to the LUCC results, the decrease in the groundwater recharge due to climate change could be even more than three times the one associated with the urban land-use growth. In any case, water security is likely threatened in the medium-term, and thee is also an increase in both flooding risk and marine intrusion [1]. In other words, the decrease in groundwater recharge as a result of LUCC is especially relevant for coastal municipalities, because of the increased associated risk of flooding and marine intrusion [27].
In the groundwater recharge evaluation, under moderate and severe climate change scenarios (RCP 4.5 and 8.5, respectively), the loss of groundwater recharge could be up to 37 times higher than the one estimated if current climate conditions continue. An increased temperature of 3 °C could generate a 28% increase in evapotranspiration. Accordingly, between 0.01% and 0.25% of the potential groundwater recharge per km 2 associated with land-use change to urban use is expected to be lost by 2050.
Furthermore, the estimated variation in potential groundwater recharge remains consistent, despite uncertainties associated with the selected hydrological model. This is only because these possible inaccuracies could be systematic and directly proportional to the groundwater baseline estimated for the initial study year.
The potential groundwater recharge was projected for different climatological conditions. We observed a greater influence of climate change than LUCC on this variation. This confirmed that both evapotranspiration increases and recharge decreases as the temperature increases due to climate change. Therefore, this work provides a methodology useful to estimating impacts on water security, mostly in low-latitude regions with high urban growth around the world.