Stochastic Climate Ensembles Reveal Hidden Uncertainty in a Depleted Mexican Aquifer

Groundwater is one of the most difficult components of the water cycle to measure directly. Unlike river discharge or rainfall, aquifer recharge cannot simply be gauged at a point; it must be inferred indirectly from rainfall, runoff, evapotranspiration, and the geology that governs how water percolates downward. In semiarid regions, where precipitation is scarce and evaporation rates are punishingly high, this inference becomes even harder. A new study published in Frontiers in Water tackles this problem head-on by building an integrated modeling framework that traces climate variability from the sky all the way down to the water table, using the overexploited Aguascalientes Valley aquifer in central Mexico as its proving ground.

The research, led by Roxana Nicte-Ha Hughes-Lomelin and colleagues at the Universidad Michoacana de San Nicolás de Hidalgo, together with a collaborator from the Institute of Geophysics at the Universidad Nacional Autónoma de México, addresses a persistent methodological gap. While scientists have long coupled surface-water and groundwater models, and separately used stochastic weather generators to simulate climate uncertainty, few frameworks explicitly propagate multiple synthetic climate realizations through a rainfall-runoff model and then into a spatially distributed groundwater-flow model. The new approach, called MASHWIN-HBV-MODFLOW, does exactly that, generating 100 equiprobable climate scenarios and following their consequences sequentially through the entire hydrological system.

The stakes in Aguascalientes could hardly be higher. The valley aquifer supplies roughly 94 percent of the state’s potable water, serving a population of nearly 1.5 million people, about 60 percent of whom live in the state capital. The region receives an average of just 510 millimeters of rain per year, while potential evaporation reaches 2,010 millimeters annually. The aquifer itself spans 3,129 square kilometers within a geological graben bounded by north-south trending normal faults, and consists mainly of sand and gravel sediments with minor silt and clay intercalations. Official figures from Mexico’s National Water Commission indicate that extraction runs at about 348 million cubic meters per year against a natural recharge of only about 250 million cubic meters, leaving a structural overdraft of nearly 100 million cubic meters annually that steadily drains the aquifer’s nonrenewable storage.

The framework’s first stage relies on MASHWIN, a multisite multivariate stochastic model developed at the Technical University of Valencia that combines autoregressive-moving-average modeling with spatial disaggregation capabilities. The researchers used it to generate synthetic monthly series of precipitation and potential evapotranspiration for seven meteorological stations, employing a periodic first-order multivariate autoregressive model. That choice reflected a deliberate statistical trade-off: for a seven-dimensional system, each additional autoregressive lag would introduce a 7-by-7 coefficient matrix, so a second-order model would require 98 lag coefficients per monthly parameterization, far too many to estimate reliably from the available 29-year observational record. Precipitation was subjected to a logarithmic transformation to reduce positive skewness, while evapotranspiration was modeled without normalization. After a five-year stochastic warm-up period was discarded, the team retained 100 synthetic realizations, each 29 years long, matching the historical record.

These synthetic climate series were then fed into the HBV rainfall-runoff model, a conceptual hydrological model that converts rain and evapotranspiration inputs into total runoff and infiltration. The HBV model was calibrated against historical runoff data from a gauged watershed and applied separately to two subbasins, Alto Aguascalientes and Medio Aguascalientes. Because independent runoff observations were insufficient for basin-specific calibration, the same parameter set was transferred unchanged to both subbasins, meaning that differences in simulated recharge between them reflect differences in meteorological forcing rather than explicit spatial representation of soil or land-use heterogeneity. The infiltration outputs, aggregated from monthly to annual values using a hydrological year running from October to September, served as the recharge input for the groundwater model.

The groundwater component used MODFLOW, configured within the ModelMuse interface as a single-layer system with 828 cells of 2 by 2 kilometers, of which 335 were active. Aquifer thickness varies from roughly 300 meters in the west to 400 and 450 meters in the east. The model ran under transient conditions over the 29-year period from 1985 to 2014, incorporating vertical recharge from HBV, groundwater extraction based on officially reported water concessions, and lateral flow across boundaries defined by the aquifer’s hydrogeological limits. Calibration against observed piezometric levels, using a split-sample strategy with 80 percent of the data for calibration and 20 percent for validation, produced strong performance metrics: a mean Nash-Sutcliffe efficiency of 0.90 during calibration and 0.83 during validation, with correlation coefficients of 0.96 and 0.92 respectively.

The results reveal a striking contrast between the two subbasins. The Alto Aguascalientes basin, despite covering approximately 2,791 square kilometers compared to Medio Aguascalientes’ 801, showed lower and more stable recharge, with a median annual recharge depth of about 2 millimeters. The smaller Medio Aguascalientes basin exhibited greater fluctuations, generally between 4 and 6 millimeters per year, with notable outliers during 2018 and 2035. When absolute infiltration volumes were normalized by basin area, the mean annual recharge depth was approximately 0.85 millimeters per year for Alto Aguascalientes versus 2.25 for Medio Aguascalientes. The researchers attribute this contrast to a combination of factors: Medio Aguascalientes is dominated by sandstone-conglomerate and rhyolitic tuff lithologies, has more extensive natural vegetation, and slightly higher proportions of medium-permeability soils, while Alto Aguascalientes features more alluvium, greater agricultural and built-up land cover, and predominantly low-permeability soils. Terrain slope, by contrast, was nearly identical between the two basins.

Perhaps the most consequential finding concerns how recharge uncertainty accumulates in the aquifer over time. Because pumping rates and lateral inflows were held identical across all 100 scenarios, the simulations isolate the effect of climate variability alone. The ensemble’s median closely tracked the deterministic simulation, with an average root-mean-square error of 21.33 meters, a Nash-Sutcliffe efficiency of 0.87, and a correlation of 0.95 between the deterministic trajectory and the ensemble’s central tendency. But the spread among individual realizations widened progressively through the simulation period. The authors interpret this as a cumulative storage response: during sequences of relatively dry years, reduced recharge deepens the annual deficit between inflow and extraction, causing additional depletion that persists across consecutive stress periods rather than resetting. Isolated wet years, meanwhile, only temporarily moderate drawdown and cannot compensate for storage losses accumulated during preceding dry spells in an aquifer already under sustained overdraft.

The authors are candid about the framework’s limitations. The single-layer groundwater model cannot resolve vertical hydraulic gradients, interlayer leakage, or confined-unconfined transitions, and the absence of explicit low-permeability interbeds may cause recharge perturbations to propagate to simulated heads more directly than they would in a multilayer system. The available dataset, they note, does not provide sufficient information to independently constrain layer geometry or layer-specific hydraulic properties, so a two-layer sensitivity model was not attempted. Piezometric observations extend only through 2014, meaning post-2014 simulations represent plausible stochastic trajectories rather than independently validated forecasts. The annual stress periods also smooth away seasonal recharge dynamics, even though roughly two-thirds of simulated infiltration in both basins occurs between May and October, and the highest-infiltration month accounts for a median of over 40 percent of the annual total.

Despite these caveats, the framework offers water managers something a single deterministic run cannot: a probabilistic baseline. Rather than relying only on the expected groundwater trajectory, managers can use the ensemble distribution to estimate the probability of low-recharge conditions, evaluate head responses under dry, median, and wet realizations, and identify periods when climatic variability raises the risk of accelerated decline. The authors suggest this could support drought contingency triggers tied to the probability of groundwater heads crossing management-defined thresholds. Looking forward, they prioritize coupling the stochastic recharge ensemble with dynamic pumping scenarios, incorporating seasonal stresses, and eventually adding multilayer hydrostratigraphy and climate-change projections. For the Aguascalientes Valley, and for the many semiarid aquifers worldwide facing similar data constraints and mounting extraction pressure, the message is clear: average conditions tell only part of the story, and the full range of plausible climate futures must be carried through the entire hydrological chain to understand what lies beneath.

Subject of Research: Stochastic propagation of climate variability through coupled surface-water and groundwater models in a semiarid Mexican aquifer

Article Title: Integrated modeling of surface water and groundwater systems: a multivariate stochastic perspective

Article References: Hughes-Lomelin, R. N.-H., Hernandez-Bedolla, J., Garcia-Perez, C. L., Hernandez-Hernandez, M. A., & Sanchez-Quispe, S. T. (2026). Integrated modeling of surface water and groundwater systems: a multivariate stochastic perspective. Frontiers in Water, 8, Article 1915442. https://doi.org/10.3389/frwa.2026.1915442

Image Credits: AI Generated

DOI: 10.3389/frwa.2026.1915442

Keywords: groundwater recharge, stochastic modeling, MASHWIN, HBV model, MODFLOW, Aguascalientes aquifer, semiarid hydrology, uncertainty propagation, rainfall-runoff modeling, Mexico, water resource management, climate variability