1. Introduction

Understanding groundwater systems in remote high mountain contexts often involves exploring and identifying aquifers within complex geological deposits (Meju, 2002), frequently with scarce to nonexistent previous information (Bishop et al., 2008). The spatial distribution and characterization of aquifer confining layers at high-altitude sites is crucial for assessing the impact of freeze-thaw cycle dynamics or long-term permafrost degradation on headwaters and consequently, on downstream ecosystems (Bense et al., 2009). For instance, changes in vertical leakage within these aquifers can lead to alterations in the mixing regimes (Ravenscroft et al., 2005) and hydraulic feedback for assessing groundwater response times, resilience, and quality of freshwater supply (Cuthbert et al., 2019; Perrone and Jasechko, 2019; Gleeson et al., 2020; Jasechko and Perrone, 2021), which are particularly critical under climate change scenarios in arid regions (Leray et al., 2013).

Assessing such issues in remote, scarcely monitored permafrost natural systems, facilitates visualizing future hydrological pathways in regions undergoing cryosphere degradation (Jorgenson et al., 2001; Jin et al., 2021). For example, high-altitude settings (>5,000 m a.s.l.) in the dry Andes (Fig. 1) enable colder, drier conditions, which under the Sixth Assessment Report (AR6) of the United Nations Intergovernmental Panel on Climate Change (IPCC) 50-year future warming scenario (IPCC, 2023), could change subsurface ice content and ground freezing duration due to snowmelt rate and partitioning (Musselman et al., 2017; Hammond et al., 2019). Consequently, groundwater flow path rerouting, affecting both subsurface confinement (Cochand et al., 2020) and annual groundwater discharge (Somers and McKenzie, 2020), is expected.

Fig. 1.  Study area in the dry Andes of northern Chile. Geophysically surveyed sites (white crosses) are shown north from the Ojos del Salado massif. The map scale decreases by about 15% from front to back. North arrow points downwards.

In frost-prone areas, such as periglacial environments, the disruption of natural flow systems upon surface warming (e.g., permafrost thaw) may impact the organization of groundwater flow components (Zinn and Konikow, 2007), as freezing temperatures modify water infiltration into the frozen subsurface (Liao and Zhuang, 2017; Uhlemann et al., 2021). Therefore, an adequate appraisal of such processes requires assessing recharge variability and active layer thicknesses, as seasonal fluctuations of unconfined groundwater levels affect aquifer volume (Bethke and Johnson, 2008) and groundwater residence times (Bethke and Johnson, 2002).

In Andean permafrost areas (Bartsch et al., 2016; Obu et al., 2019), such as the Ojos del Salado massif (OSM) in northern Chile (27°06’ S, 68°32’ W; Fig. 1), future thaw of interstitial ice could impact aquifer recharge, as near surface, low-permeability permafrost conditions strongly affect the distribution of hydraulic heads and associated fluid flow patterns (e.g., Bense and Person, 2008). Therefore, the role of permafrost in mountain hydrology should ideally be more relevant to water budget (Arenson et al., 2022), as permafrost thaw results in the deepening of aquitard roofs, allowing deeper flow paths and causing a decrease in summer stream temperatures and evapotranspiration in headwater catchments (Sjöberg et al., 2021).

The present work aims to unveil permafrost structures in high-altitude aquifers along an altitudinal gradient (4,550-5,830 m a.s.l.) at the OSM, a cryosphere-relevant massif in the dry Andes of South America, which may present different cryo-hydrogeological configurations (Fig. 2) and other hillslope permafrost environments (Evans and Ge, 2017). Different locations were assessed along a NNW-SSE transect, from an inferred no-permafrost condition at the warmest, lowest altitude site (4,550 m a.s.l.), an inferred sporadic permafrost condition at an intermediate altitude site (5,260 m a.s.l.), and a positive control for permafrost at the coldest, highest site (5,830 m a.s.l.) (Fig. 2). A standard methodology was used to characterize permafrost areas by contrasting electrical resistivities outlining the active layer depth (seasonal ground thaw) and permafrost table (perennially frozen ground) to interpret the current relevant subsurface architecture.

Fig. 2. Conceptual subsurface structures at the Ojos del Salado transect (see Fig. 1 for location of sites). Arrows represent idealized flow paths. A. Highest-altitude, colder configuration with consistent, homogeneous permafrost distribution. B. Configuration with isolated permafrost distribution and possibly the lower permafrost limit. C. Warmer areas, water leakage and negative permafrost control configuration.

2. Study area and methods

The study area was chosen since there are neither permafrost altitudinal zonation assessments nor subsurface characterization of permafrost thicknesses. At each site, electrical resistivity tomography (ERT) data (Krautblatter and Hauck, 2007; Hilbich et al., 2021; Buckel et al., 2022) were integrated along with permafrost probability estimates (Ran et al., 2015) derived from long-term ground temperature datasets (Fig. 3). Such a geophysical survey covering different altitudes at sites with decade-long surface temperature data is unprecedented in the Andes.

Fig. 3. Geophysically surveyed areas at the Ojos del Salado massif (A: Murray; B: Atacama; C: Tejos). Electrical resistivity tomography transects (green circles, white lines) and water sample locations (blue triangles) are shown for each site. The Atacama site has two perpendicular transects: a longer (195 m) N-S and a shorter (95 m) E-W. Contour lines every 25 m.

2.1. Characteristics of the study area

The study area is located around the Ojos del Salado volcano (6,893 m a.s.l.), which is the highest massif of the dry Andes. The Ojos del Salado massif (OSM) and its immediate surroundings are part of the mountain desert and tundra belt of the high Andes of northern Chile. Therefore, the massif lacks active glaciers, showing only perennial snow/firn patches. The thermal climatic snowline runs at an altitude of about 7,000 m a.s.l. (Houston and Hart, 2004), and precipitation is dominated by thin winter snow cover only for a few weeks between May and October (Kereszturi et al., 2022). The massif has a predicted permafrost condition (>0.5 probability) above 5,000 m a.s.l. (Bartsch et al., 2016). Additionally, thermokarst features near the Atacama site at 5,260 m a.s.l. are present, including buried ice layers (Kereszturi et al., 2022).

Regolith thicknesses around the established transects account for several meters of coarse-grained volcanic debris. Wind erosion, transport, and deposition are the most active surface processes at the study sites (Nagy et al., 2019), where the boulders show marked features of wind abrasion. Evidence of periglacial cryoturbation and frost heave is absent on these surfaces, while slope processes are present only above ~5,600 m a.s.l. (Nagy et al., 2019).

According to temperature monitoring data between 2012 and 2023, the study sites show an overall warming soil temperature trend (Fig. 4). In fact, the mean warming was about 0.03 °C/month at the Murray site, 0.012 °C/month at the Atacama site, and 0.006 °C/month at the Tejos site. Satellite-based snow cover studies show no winter cover at Murray and a predominant winter cover at both Atacama and Tejos, with the latter presenting an occasional snow cover in mid-summer (Nagy et al., 2019, 2023).

The composition and characteristics of the substrate for each study site are explained below. We provide bulk density and porosity estimates for the three sites and relative and absolute humidity values for Murray.

Fig. 4. Daily temperature time series at Murray (blue line), Atacama (red line), and Tejos (yellow line) for three different depths, during a 9- to 11-year period. Ground surface temperature datasets extracted from the PermaChile network (www.permachile.com).

2.1.1. Murray site

The Murray site is located ~20 km northwards from the OSM at 4,550 m a.s.l., in a pediment surface developed on the central part of a fossil alluvial fan. There are neither watercourses nor evident drainage lines around the site. The largest clasts on the surface are cm-sized pumice fragments. The bulk density and porosity of the regolith were determined following Nagy et al. (2019), and resulted in 1.4 g/cm³ at 10 cm, 1.5 g/cm³ at 35 cm, and 1.4 g/cm³ at 60 cm depth, and 21-25% v/v for the same depth range (the lowest porosities of all sites surveyed). The relative humidity was measured by means of a HOBO Pro v2 Temperature/Relative Humidity logger on an hourly basis. At 10 cm depth, it drops to 47% in summer with an absolute humidity value of 3.75 g/m³.

2.1.2. Atacama site

The Atacama site is located ~5.5 km northwards from the OSM at 5,260 m a.s.l., in a flat lava plain. It is only above this site that slopes become relevant (see Fig. 1). The site is characterized by extensive megaripple fields, which extend across the rocky outcrop interspersed with large, wind-carved boulders (Nagy et al., 2019, 2020). The surveyed area is rich in pumice fragments with a well-developed lag gravel pavement, and at about 35 cm depth the regolith’s >2 mm diameter fraction is dominant. In the 60 cm depth range, the porosity of the sediment is 50-55% v/v and the water absorption capacity (calculated following Nagy et al., 2019) is high only in the upper 10 cm (1,081 nm) and then drastically decreases with depth. The bulk density of the regolith is 1.3 g/cm³ at a depth of 10 cm and 1.5 g/cm³ at depths of 35 and 60 cm. Below 40 cm depth, rock fragments larger than 10 cm appear and the sediment becomes highly compact. Evidence of cryoturbation was not observed (Nagy et al., 2019).

At this site, the phase-change model of Nagy et al. (2020) applied to the near-surface (first 100 cm depth) regolith showed that the active layer is at least >1 m-thick. No evident ice was detected, but the model suggested that permafrost was a possibility. Mean annual temperatures are typically close to 0 °C, although warmer temperatures were sometimes recorded at all depths (Table 1). The number of days with temperatures above 0 °C decreases with depth and increases at lower altitudes (Table 1). These temperature data indicate that the presumed subsurface permafrost at Atacama may play a role in keeping the regolith frozen during specific periods. At a depth of 60 cm, the proportion of freeze-thaw cycle days (i.e., days when the air temperature fluctuates between freezing and non-freezing temperatures) estimated at Atacama decreases to a third of that estimated at Murray (Table 1).

According to Nagy et al. (2019, 2020), thawing starts in mid-October and the surface freezes again in early to mid-April. It takes 40-50 days to thaw the first 100 cm and refreezing is delayed by about the same amount of time. Typically, the subsurface does not fully refreeze until late May or the first half of June. The onset of the thawing period may be delayed for up to one month in case of persisting snow cover, although the end of the thawing period rarely shifts as much, as by late summer there is no significant snowfall or snow accumulation.

2.1.3. Tejos site

The Tejos site is located ~2.5 km northwards from the OSM at 5,830 m a.s.l., in a mountain tundra periglacial environment. The site is in a flat area amid thick lava flows, moraine ridges, and rock glacier deposits.

The site is part of a perennial firn-fed outwash plain, with occasional snow and wind-carved boulders without sand accumulations or traces of cryoturbation. The slopes are dominated by solifluction lobes, with huge debris slopes leading to the top of the massif. The specific monitoring location does not show signs of meltwater flooding but a consistent blanket of gravel lag deposits, with the largest grain size out of the three monitoring sites. As at Atacama, the regolith’s coarser (>2 mm) fraction is finer at 35 cm depth, with a higher silt fraction around this depth as well. The bulk density of the regolith is the lowest out of all sites: 1.1, 1.3, and 1.0 g/cm³ at depths of 10, 35, and 60 cm, respectively. The porosity of the sediment is even higher than at Atacama, ranging between 55-67% v/v (Nagy et al., 2019).

At this site, the daily temperature measurements (Table 1) show a persistent presence of ice-bearing permafrost at least since February 2012. Thermal loggers placed near the permafrost table at 60 cm depth are occasionally frozen until mid-February. According to the phase-change model of Nagy et al. (2020), the active layer thickens to a maximum of 70-80 cm. At the surface, regolith thawing starts in early November and freezes again in early April. The thawing reaches the bottom of the active layer by the end of February or the beginning of March, i.e., requiring between 100-120 days. The refreezing rate at Tejos is twice as fast as at Atacama. The onset of the thawing period can be delayed by up to a month due to the snow cover effect and its termination shortened to the beginning of March, although transient, refreezing events may occur at any stage of the melting period particularly in the upper 10-20 cm. This delayed thawing occurs when snow does not sublimate following summer snowfalls but rather it melts over several weeks, with meltwater soaking the surface, creating temporary ice cement.

2.2. Determination of subsurface permafrost structures

2.2.1. Freezing and thawing indicators

The likelihood of permafrost occurrence in the three study sites was analyzed by means of surface temperature datasets from 2012-2023 for the Atacama and Tejos sites, and 2014-2023 for the Murray site (Fig. 4 and Table 1). All thermal loggers were located on a horizontal surface, far from any shading topographic feature (see Nagy et al., 2019 for more information on the measurement procedure). The loggers were of the type HOBO Pro v2 (U22-001), with an operation range from -40 to 70 °C, an accuracy of ±0.21 °C from 0 to 50 °C, and a resolution of 0.02 °C at 25 °C.

The thermal exposure was estimated for each site through the freezing/thawing degree days and the Frost number (F⁺) (Nelson and Outcalt, 1987; Barry and Gan, 2011, p. 172). The total annual freezing (FDD) and thawing (TDD) days are defined as the cumulative number of days in a year with daily air temperatures below and above 0 °C, respectively. The Frost number is a simplified index for the likelihood of permafrost occurrence and is calculated as follows:

Eq. 1

These estimates compare the overall magnitude of thermal exposure between sites, such as the n-factor determination (see Klene et al., 2001), and help compare frost exposure and permafrost conditions in periglacial environments (Ran et al., 2015).

Frost numbers can define the occurrence probability of permafrost as follows: continuous (F⁺≥0.67), discontinuous (0.67>F⁺≥0.6), sporadic (0.6>F⁺ ≥0.5), and no (F⁺<0.5) permafrost (Barry and Gan, 2011).

2.2.2. Electrical resistivity tomography

Electrical resistivity tomography (ERT) is a standard geophysical method to investigate subsurface structures and identify their properties. By combining ERT with permafrost occurrence estimates, more confident permafrost thicknesses can be determined (Nagy et al., 2020).

All ERT surveys conducted during this study took place around noon. The data were collected with a 4-point light Lippmann resistivity meter with a distance between electrodes of 10 m for the Murray and Atacama sites and 5 m for the Tejos site. The array employed for the resistivity lecture was dipole-dipole (see Kneisel, 2006) due to a better performance in terms of penetration depth. Four profiles were surveyed, one at Murray (195 m length), two at Atacama (95 and 195 m length), and one at Tejos (145 m length) (Table 2). To process the data and generate the calculated resistivity profiles, the Res2DInv 5.0 inversion software (Loke and Barker, 1996) was used. In all ERT surveys, there were no significant error issues on electrode contact resistance.

3. Results

3.1. Permafrost occurrence

The Frost number estimations were 0.34, 0.53, and 0.72 for the Murray, Atacama, and Tejos sites, respectively. Surface temperature datasets at 10 cm depth evidence frost conditions representative of a permafrost environment at Tejos, a sporadic occurrence at Atacama and no permafrost at Murray. These estimates suggest an altitudinal gradient for mountain permafrost occurrence at the OSM, accounting for a negative control, a fringe area, and a favorable site for frozen grounds by the end of the warm season.

3.2. Geophysical data

Figure 5 presents the electrical resistivity profiles obtained at the three sites. The ERT transect at Murray (Fig. 5A) shows two layers and a half-space separated by two evident resistivity contrasts. First, a ~2-6 m-thick surface layer with low resistivity values (0.2-0.5 kΩm), beneath it, a meter-thick, intermediate layer not parallel to the surface with higher resistivity values (1-7 kΩm) and a half-space with resistivities <5 kΩm below.

Fig. 5. Electrical resistivity profiles. A. Murray, N-S transect. B. Atacama, N-S transect. C. Atacama, W-E transect. D. Tejos, N-S transect. See figure 3 for a location map. Black triangles are 5 m scale for reference, not the electrode spacing. Isocontours are given in Ωm.

At Atacama, the N-S transect (Fig. 5B) evidences a ~3 m-thick surface layer with resistivities ranging between 0.5-2 kΩm. A clear high resistivity layer (15-20 kΩm) appears in the southern half of the profile and thins northwards. In the E-W resistivity profile (Fig. 5C), it is possible to define three different resistivity layers: (i) a shallow, 3 m-thick, low resistivity layer (0.5-2.5 kΩm); (ii) an intermediate, 20 m-thick, high resistivity layer (15-20 kΩm); and (iii) a bottom half-space with resistivities lower than 5 kΩm.

At Tejos, there is a shallow, 2-5 m-thick layer with resistivities below 10 kΩm, that overlies a 20 m-thick domain with resistivity values of up to 100 kΩm. The Tejos site exhibits higher resistivity domains than those observed at the other sites and more contrasting resistivity transitions overall.

4. Discussion

The thermal analysis evidences a transition from mountain desert without permafrost (Murray) to a mountain periglacial zone with ice-rich permafrost (Tejos), where the active layer is assumed to be present above ~5,260 m a.s.l. In fact, at Tejos, the thawing degree days are dramatically lower when compared to Atacama and Murray, even though the percentage of freeze-thaw cycle days does not show this contrast as evidently. This aspect is clearly shown in Figure 4 as well, where temperatures at 60 cm depth remained always below 0 °C at Tejos, indicative of either a permafrost-table top or a thermal transition for ice-phase change.

4.1. Discretization of permafrost structures

The ERT results do not unequivocally confirm an altitudinal gradient from no frozen grounds (Murray) up to subsurface values common for permafrost bodies (Tejos). However, there is a consistent resistivity increase at Atacama below 3 m depth, reaching values between 15-20 kΩm. Below 20 m depth, there are important differences in the resistivity values and their spatial extent for all sites. For example, at Murray, there is a large continuous body <5 kΩm, while at Atacama there are lens-like bodies with resistivities of up to 20 kΩm above a less homogeneous <5 kΩm domain. The high-resistive features in the Atacama site thin northwards and longitudinally (Fig. 5C), while at Tejos, in contrast, the predominant subsurface feature is a highly resistive body reaching even 100 kΩm in certain parts.

According to the geological information of the study area (Clavero et al., 2012; Naranjo et al., 2019), the high resistivity values measured at Murray can be indicative of the Laguna Verdeignimbrite. At shallow depths, the low resistivity values could refer to salty alluvial sediments whereas the high-resistivity m-thick body beneath could be a confined aquifer, possibly part of a sedimentary layer just above the ignimbrite. At higher altitudes, ERT data at the Atacama and Tejos sites show the presence of >10 kΩm resistivity structures in the subsurface, possibly lava flows, whose spatial variability can be explained by the phase change of the interstitial water. The ~10 kΩm threshold can particularly be related to the limit between frozen and no frozen material. At the Tejos site, the <5 m-thick, low resistivity shallow body would be indicative of alluvial deposits.

The high resistivity domains imaged particularly at Tejos are consistent with those observed at around 30-31° S in the Argentine Andes (Halla et al., 2021; Villarroel et al., 2022) and at around 28° S in the Chilean Andes (Hilbich et al., 2021; Mathys et al., 2022). These high resistivity values are also consistent with findings in the Alps (Krautblatter and Hauck, 2007; Buckel et al., 2022), establishing unfrozen states below ~10 kΩm and comparable to measured rock resistivity values near 0 ºC (Scandroglio et al., 2021). The bodies imaged in the Alps may be associated with a putative transition layer between the active layer and the permafrost table (Arenson et al., 2022). Therefore, the Tejos high resistivity domains are interpreted as frozen materials with an indeterminate quantity of ice.

By taking together ERT analyses, subsurface frost depths, and Frost numbers, we suggest that the most probable cause for the order-of-magnitude increase in electrical resistivity from Murray to Atacama and Tejos is the increasing presence of frozen ground.   It is also worth noticing that as the ERT survey took place at the end of summer, the high resistivity anomalies detected cannot be snowstorm driven.

4.2. Cryo-hydrogeological implications

As the active layer stores and conducts water (Fig. 2), the thawing period and depth are indicators of a warming environment impacting surface hydrology (Schuur et al., 2015; Walvoord and Kurylyk, 2016). For example, a deeper thawing implies a deepening of flow paths, expanding unconfined aquifers (Sjöberg et al., 2021) and perturbing the background flow field (Wang et al., 2014), therefore enabling transitions into different hydrogeological frameworks (Fig. 2).

The assumption of a continuous transition layer at around 10 kΩm, limiting frozen from unfrozen sediments, implies a consistent frozen body from the shallow surface down to at least 25 m depth as in the case of the Tejos site (Fig. 5D). At Atacama, the same resistivity threshold criterion suggests a lenticular and wedge-like shape for these frozen bodies, with thicknesses below 10 m (Fig. 5B, C). The presence of a consistent, rather homogeneous active layer at Tejos implies that upon surface recharge, the sediment transport yield could accumulate finer grain size above the transition layer. In contrast, the lens-like distribution observed at Atacama could imply a less effective sediment transport with multi-entry points. The potential finding of frozen, lens-like features with distinct salinity (as in deserts) and grain size could represent moderately conductive anomalies along resistive profiles, or a different case with partially frozen sediments (Campbell et al., 2021) with resistivities portraying the ice-phase change transition.

In arid areas, high surface salinity results from former water bodies affected by cryogenic desiccation and wind activity (Kereszturi et al., 2022). A ground with high salinity will have a low resistivity when imaged by ERT techniques, which suggests that some of the low resistivities observed in the study area would be indicative of high saline conditions in the shallow ground. Ground surface salinity measured at Atacama (78.1 µS/cm) and Tejos (224 µS/cm) could represent the effect of in-situ salinity concentration (Aszalós et al., 2020).

In volcanic areas, preferential pathways for fluid circulation exist due to a high hydraulic conductivity gradient between structures (Vittecoq et al., 2019). In addition, geothermal gradient studies in volcanic areas require basic thermal modeling to contrast each site’s structural and geophysical data (Nagy et al., 2020). Hence, despite considering a high geothermal gradient as the bottom boundary condition, a bottom-up thawing cannot be longer sustained based on the ERT data, as results are highly indicative of the presence of m-thick frozen structures in the subsurface.

4.3. Limitations and research opportunities

Geophysical surveys around poorly monitored areas without boreholes or previous geophysical assessments are challenging, particularly in high-altitude, remote areas. As such, uncertainties emerge and need to be identified and quantified accordingly. The use of analog hydrogeological settings can be a useful approach. Analogs to the case study presented here can be found in Martian interpretations (Michalski et al., 2013), where volcanically outgassed water would have been locked as taliks, the water equivalent being periodically mobilized, recharging the subsurface through basal thawing.

One limitation of the ERT method is the lack of statistical calibration for the inferred subsurface materials. This deals with both the pore-size-related complexity in saturation estimation and the reactance of either liquid water or ice, as both materials differ in their dielectric properties. For instance, when working over resistive material, the interpreted electrical resistance should account for the ubiquitous presence of weakly conductive and nonconductive pores in rocks and sediments (Zhu et al., 2023). In permafrost areas, these bodies could behave either as Archie or non-Archie rocks depending on their clay content and temperature. Furthermore, in the case of interstitial water/ice, concentration and uniformity determinations need to be calibrated if water equivalent data is required.

Another source of uncertainty is the possible phase transitions from thawed to frozen ground in fine-grained and saline sediments (O’Neill et al., 2019). This can be addressed by identifying the resistivity values along the transition layers that surround highly resistive bodies, allowing more precise geometrical constraints and water-equivalent quantification. By taking this approach, it may be possible to track permafrost degradation in terms of pore water equivalent in larger areas, therefore providing more applied insights for upscaling the effects of climate change on mountainous environments.

Groundwater flow modeling in permafrost settings often lacks robust hydrogeological field data (Kurylyk and Walvoord, 2021, p. 514), so the assessment presented here aims to strengthen the physical foundations for high Andean environments by providing geocryological data for subsequent evaluations of critical groundwater thresholds. Lastly, from a bottom-up perspective, potential recharge impacts on vadose zone infiltration remained unaccounted for. It is therefore essential to consider that, when dealing with volcanic settings, the computation of geothermal gradients should consider: (i) temperature boreholes deeper than 20 m; (ii) an indirect characterization of temperatures from ERT profiles against resistivity curves of thawing sediments; and (iii) hydrochemical analyses as a proxy of water subsurface origin with temperature variability as a proxy of depth. Overall, direct temperature and/or geotechnical data are required to constrain the inversion models and reduce the uncertainty in the derived results.

5. Conclusions

This study confirms a permafrost altitudinal gradient between 4,550-5,830 m a.s.l. at the Ojos del Salado massif in the dry Andes by integrating electrical resistivity surveys and decade-long ground surface temperature datasets.

At sites where the temperature data favors permafrost occurrence, the resistivity surveys reveal consistent differences attesting to high resistivities associated with permafrost thicknesses of about 25 m in the coldest and highest altitude site (Tejos). Abrupt resistivity gradients and well-defined ~10 kΩm resistivity bodies at Atacama and Tejos are indicative of possible transition layers. These results are consistent with previous studies on active rock glaciers in the Alps and the Andes and confirm the presence of frozen ground layers at the Ojos del Salado massif.

This study has implications for better understanding groundwater recharge, surface water dynamics, and the overall response of hydrological systems to climate change in high-altitude, cold-dry settings. Despite their limitations, geophysical surveys afford research opportunities for cold, remote sites. Further integration of field data and improved modeling approaches will be crucial for evaluating groundwater dynamics, groundwater-surface water partition and their effects on streamflow, as well as future tipping points of stagnation, upwelling events, flow path variation, and novel hydrological connections in these sensitive environments.

Acknowledgments
Author contributions. SRP and SL organized the first manuscript drafts, fieldwork, and data analysis. GY contributed to geophysical formal analysis and interpretation. SRP, SL, and FS analyzed and interpreted thermal data and worked on the final paper focus, structure, and organization. SRP, EM, and EB performed fieldwork. BN integrated field temperature data and participated in fieldwork at the Atacama and Tejos sites. This work was supported by Beca Postdoctorado Escuela de Ingeniería, Pontificia Universidad Católica de Chile, and the ANID-FONDECYT Postdoctorado Nacional Folio 3230723. The authors also acknowledge the Laboratorio de Evaluación y Control de Calidad del Agua, DIHA-PUC, Chile. All authors reviewed and actively commented on the manuscript. Alfonso Fernández and the Editor helped reviewing this manuscript.

References

Arenson, L.U.; Harrington, J.S.; Koenig, C.E.M.; Wainstein, P.A. 2022. Mountain permafrost hydrology-a practical review following studies from the Andes. Geosciences 12 (2): 48. https://doi.org/10.3390/geosciences12020048

Aszalós, J.M.; Szabó, A.; Megyes, M.; Anda, D.; Nagy, B.; Borsodi, A.K. 2020. Bacterial diversity of a high-altitude permafrost thaw pond located on Ojos del Salado (Dry Andes, Altiplano-Atacama Region). Astrobiology 20 (6): 754-765. https://doi.org/10.1089/ast.2018.2012

Barry, R.; Gan, T.Y. 2011. The global cryosphere: past, present and future. Cambridge University Press: 472 p. https://doi.org/10.1017/CBO9780511977947

Bartsch, A.; Grosse, G.; Kääb, A.; Westermann, S.; Strozzi, T.; Wiesmann, A.; Duguay, C.; Seifert, F.M.; Obu, J.; Goler, R. 2016. GlobPermafrost-how space-based Earth observation supports understanding of permafrost. In Living Planet Symposium, Proceedings, ESA SP: 740 p. Prague.

Bense, V.F.; Ferguson, G.; Kooi, H. 2009. Evolution of shallow groundwater flow systems in areas of degrading permafrost. Geophysical Research Letters 36 (22): L22401. https://doi.org/10.1029/2009GL039225

Bense, V.F.; Person, M.A. 2008. Transient hydrodynamics within intercratonic sedimentary basins during glacial cycles. Journal of Geophysical Research, Earth Surface 113 (F4): F04005. https://doi.org/10.1029/2007JF000969

Bethke, C.M.; Johnson, T.M. 2002. Paradox of groundwater age. Geology 30 (2): 107-110. https://doi.org/10.1130/0091-7613(2002)030<0107:POGA>2.0.CO;2

Bethke, C.M.; Johnson, T.M. 2008. Groundwater age and groundwater age dating. Annual Review of Earth and Planetary Sciences 36: 121-152. https://doi.org/10.1146/annurev.earth.36.031207.124210

Bishop, K.; Buffam, I.; Erlandsson, M.; Folster, J.; Laudon, H.; Seibert, J.; Temnerud, J. 2008. Aqua Incognita: the unknown headwaters. Hydrological Processes 22 (8): 1239-1242. https://doi.org/10.1002/hyp.7049

Buckel, J.; Mudler, J.; Gardeweg, R.; Hauck, C.; Hilbich, C.; Frauenfelder, R.; Kneisel, C.; Buchelt, S.; Blöthe, J.H.; Hördt, A. 2022. Identifying mountain permafrost degradation by repeating historical electrical resistivity tomography (ERT) measurements. The Cryosphere 17 (7): 2919-2940. https://doi.org/10.5194/tc-17-2919-2023

Campbell, S.W.; Briggs, M.; Roy, S.G.; Douglas, T.A.; Saari, S. 2021. Ground-penetrating radar, electromagnetic induction, terrain, and vegetation observations coupled with machine learning to map permafrost distribution at Twelvemile Lake, Alaska. Permafrost and Periglacial Processes 32 (3): 407-426. https://doi.org/10.1002/ppp.2100

Clavero, J.; Mpodozis, C.; Gardeweg, M.; Valenzuela, M. 2012. Geología de las áreas Laguna Wheelwright y Paso San Francisco, Región de Atacama. Servicio Nacional de Geología y Minería, Carta Geológica de Chile, Serie Geología Básica 139-140: 32 p. Santiago.

Cochand, M.; Molson, J.; Barth, J.A.C.; van Geldern, R.; Lemieux, J.M.; Fortier, R.; Therrien, R. 2020. Rapid groundwater recharge dynamics determined from hydrogeochemical and isotope data in a small permafrost watershed near Umiujaq (Nunavik, Canada). Hydrogeology Journal 28: 853-858. https://doi.org/10.1007/s10040-020-02109-x

Cuthbert, M.O.; Gleeson, T.; Moosdorf, N.; Befus, K.M.; Schneider, A.; Hartmann, J.; Lehner, B. 2019. Global patterns and dynamics of climate-groundwater interactions. Nature Climate Change 9 (2): 137-141.

Evans, S.G.; Ge, S. 2017. Contrasting hydrogeologic responses to warming in permafrost and seasonally frozen ground hillslopes. Geophysical Research Letters 44 (4): 1803-1813. https://doi.org/10.1002/2016GL072009

Gleeson, T.; Cuthbert, M.; Ferguson, G.; Perrone, D. 2020. Global groundwater sustainability, resources, and systems in the Anthropocene. Annual review of earth and planetary sciences 48 (1): 431-463.

Halla, C.; Blöthe, J.H.; Tapia Baldis, C.; Trombotto Liaudat, D.; Hilbich, C.; Hauck, C.; Schrott, L. 2021. Ice content and interannual water storage changes of an active rock glacier in the dry Andes of Argentina. The Cryosphere 15 (2): 1187-1213. https://doi.org/10.5194/tc-15-1187-2021

Hammond, J.C.; Harpold, A.A.; Weiss, S.; Kampf, S.K. 2019. Partitioning snowmelt and rainfall in the critical zone: effects of climate type and soil properties. Hydrology and Earth System Sciences 23 (9): 3553-3570. https://doi.org/10.5194/hess-23-3553-2019

Hilbich, C.; Hauck, C.; Mollaret, C.; Wainstein, P.; Arenson, L.U. 2021. Towards accurate quantification of ice content in permafrost of the Central Andes, Part 1: geophysics-based estimates from three different regions. The Cryosphere 16 (5): 1845-1872. https://doi.org/10.5194/tc-16-1845-2022

Houston, J.; Hart, D. 2004. Theoretical head decay in closed basin aquifers: an insight into fossil groundwater and recharge events in the Andes of northern Chile. Quarterly Journal of Engineering Geology and Hydrogeology 37 (2): 131-139. https://doi.org/10.1144/1470-9236/04-007

Jasechko, S.; Perrone, D. 2021. Global groundwater wells at risk of running dry. Science 372 (6540): 418-421. https://doi.org/10.1126/science.abc2755

Jin, X.-Y.; Jin, H.-J.; Iwahana, G.; Marchenko, S.S.; Luo, D.-L.; Li, X.-Y.; Liang, S.-H. 2021. Impacts of climate-induced permafrost degradation on vegetation: A review. Advances in Climate Change Research 12 (1): 29-47. https://doi.org/10.1016/j.accre.2020.07.002

Jorgenson, M.T.; Racine, C.H.; Walters, J.C.; Osterkamp, T.E. 2001. Permafrost degradation and ecological changes associated with a warming climate in central Alaska. Climatic Change 48: 551-579. https://doi.org/10.1023/A:1005667424292

Kereszturi, A.; Aszalós, J.; Heiling, Zs.; Ignéczi, A.; Kapui, Zs.; Király, Cs.; Leél-Őssy, Sz.; Szalai, Z.; Nemerkényi, Zs.; Pál; B.; Skultéti, A.; Nagy, B. 2022. Wind-snow interactions at the Ojos del Salado region as a potential Mars analogue site in the Altiplano-Atacama desert region. Icarus 378. https://doi.org/10.1016/j.icarus.2022.114941

Klene, A.E.; Nelson, F.E.; Shiklomanov, N.I.; Hinkel, K.M. 2001. The n-factor in natural landscapes: variability of air and soil-surface temperatures, Kuparuk River Basin, Alaska, USA. Arctic, Antarctic, and Alpine Research 33 (2): 140-148. https://doi.org/10.1080/15230430.2001.12003416

Kneisel, C. 2006. Assessment of subsurface lithology in mountain environments using 2D resistivity imaging. Geomorphology 80 (1-2): 32-44. https://doi.org/10.1016/j.geomorph.2005.09.012

Krautblatter, M.; Hauck, C. 2007. Electrical resistivity tomography monitoring of permafrost in solid rock walls. Journal of Geophysical Research, Earth Surface 112 (F2): F02S20. https://doi.org/10.1029/2006JF000546

Kurylyk, B.L.; Walvoord, M.A. 2021. Permafrost hydrogeology. In Arctic hydrology, permafrost and ecosystems (Yang, D.; Kane, D.L.; editors). Springer: 493-523. Cham. https://doi.org/10.1007/978-3-030-50930-9_17

Leray, S.; de Dreuzy, J.-R.; Bour, O.; Bresciani, E. 2013. Numerical modeling of the productivity of vertical to shallowly dipping fractured zones in crystalline rocks. Journal of Hydrology 481: 64-75. https://doi.org/10.1016/j.jhydrol.2012.12.014

Liao, C.; Zhuang, Q. 2017. Quantifying the role of permafrost distribution in groundwater and surface water interactions using a three-dimensional hydrological model. Arctic, Antarctic, and Alpine Research 49 (1): 81-100. https://doi.org/10.1657/AAAR0016-022

Loke, M.H.; Barker, R.D. 1996. Rapid least-squares inversion of apparent resistivity pseudosections by a quasi-Newton method. Geophysical Prospecting 44 (1): 131-152. https://doi.org/10.1111/j.1365-2478.1996.tb00142.x

IPCC (Intergovernmental Panel on Climate Change) 2023. Summary for Policymakers. In Climate Change 2023: Synthesis Report. Contribution of Working Groups I, II and III to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (Core Writing Team; Lee, H.; Romero,J.; editor). Intergovernmental Panel on Climate Change: 1-34. Geneva. https://doi.org/10.59327/IPCC/AR6-9789291691647.001

Mathys, T.; Hilbich, C.; Arenson, L.U.; Wainstein, P.A.; Hauck, C. 2022. Towards accurate quantification of ice content in permafrost of the Central Andes-Part 2: an upscaling strategy of geophysical measurements to the catchment scale at two study sites. The Cryosphere 16 (6): 2595-2615. https://doi.org/10.5194/tc-16-2595-2022

Meju, M.A. 2002. Geoelectromagnetic exploration for natural resources: models, case studies and challenges. Surveys in Geophysics 23: 133-206. https://doi.org/10.1023/A:1015052419222

Michalski, J.R.; Cuadros, J.; Niles, P.B.; Parnell, J.; Deanne Rogers, A.; Wright, S.P. 2013. Groundwater activity on Mars and implications for a deep biosphere. Nature Geoscience 6: 133-138. https://doi.org/10.1038/ngeo1706

Musselman, K.N.; Clark, M.P.; Liu, C.; Ikeda, K.; Rasmussen, R. 2017. Slower snowmelt in a warmer world. Nature Climate Change 7: 214-219. https://doi.org/10.1038/nclimate3225

Nagy, B.; Ignéczi, Á.; Kovács, J.; Szalai, Z.; Mari, L. 2019. Shallow ground temperature measurements on the highest volcano on Earth, Mt. Ojos del Salado, Arid Andes, Chile. Permafrost and Periglacial Processes 30(1): 3-18. https://doi.org/10.1002/ppp.1989

Nagy, B.; Kovács, J.; Ignéczi, Á.; Beleznai, S.; Mari, L.; Kereszturi, Á.; Szalai, Z. 2020. The thermal behavior of ice-bearing ground: the highest cold, dry desert on Earth as an analog for conditions on Mars, at Ojos del Salado, Puna de Atacama-Altiplano region. Astrobiology 20 (6): 701-722. https://doi.org/10.1089/ast.2018.2021

Nagy, B.; Ignéczi, Á.; Kovács-Székely, I.; Pereira, S.R.; Mihajlik, G.; Felkai, P.; Mari, L. 2023. The challenges of commercial mountaineering on the highest Volcanic Seven Summit, the Ojos del Salado. Hungarian Geographical Bulletin 72 (1): 23-40. https://doi.org/10.15201/hungeobull.72.1.2

Naranjo, J.A.; Hevia, F.; Arcos, R.; Polanco, E. 2019. Geología de las áreas Nevado Ojos del Salado y Cerro El Fraile, Región de Atacama. Servicio Nacional de Geología y Minería, Carta Geológica de Chile, Serie Geología Básica 147: 37 p. Santiago.

Nelson, F.E.; Outcalt, S.I. 1987. A computational method for prediction and regionalization of permafrost. Arctic and Alpine Research 19 (3): 279-288. https://doi.org/10.1080/00040851.1987.12002602

Obu, J.; Westermann, S.; Kääb, A.; Bartsch, A. 2019. Ground Temperature Map, 2000-2016, Andes, New Zealand and East African Plateau Permafrost (dataset). University of Oslo, PANGAEA. https://doi.org/10.1594/PANGAEA.905512

O’Neill, H.B.; Wolfe, S.A.; Duchesne, C. 2019. New ground ice maps for Canada using a paleogeographic modelling approach. The Cryosphere 13 (3): 753-773. https://doi.org/10.5194/tc-13-753-2019

Perrone, D.; Jasechko, S. 2019. Deeper well drilling an unsustainable stopgap to groundwater depletion. Nature Sustainability 2: 773-782. https://doi.org/10.1038/s41893-019-0325-z

Ran, Y.; Li, X.; Jin, R.; Guo, J. 2015. Remote sensing of the mean annual surface temperature and surface frost number for mapping permafrost in China. Arctic, Antarctic, and Alpine Research 47 (2): 255-265. https://doi.org/10.1657/AAAR00C-13-306

Ravenscroft, P.; Burgess, W.G.; Ahmed, K.M.; Burren, M.; Perrin, J. 2005. Arsenic in groundwater of the Bengal Basin, Bangladesh: distribution, field relations, and hydrogeological setting. Hydrogeology Journal 13: 727-751. https://doi.org/10.1007/s10040-003-0314-0

Scandroglio, R.; Draebing, D.; Offer, M.; Krautblatter, M. 2021. 4D quantification of alpine permafrost degradation in steep rock walls using a laboratory‐calibrated electrical resistivity tomography approach. Near Surface Geophysics 19 (2): 241-260. https://doi.org/10.1002/nsg.12149

Schuur, E.A.G.; McGuire, A.D.; Schädel, C.; Grosse, G.; Harden, J.W.; Hayes, D.J.; Hugelius, G.; Koven, C.D.; Kuhry, P.; Lawrence, D.M.; Natali, S.M.; Olefeldt, D.; Romanovsky, V.E.; Schaefer, K.; Turetsky, M.R.; Treat, C.C.; Vonk, J.E. 2015. Climate change and the permafrost carbon feedback. Nature 520: 171-179. https://doi.org/10.1038/nature14338

Sjöberg, Y.; Jan, A.; Painter, S.L.; Coon, E.T.; Carey, M.P.; O’Donnell, J.A.; Koch, J.C. 2021. Permafrost promotes shallow groundwater flow and warmer headwater streams. Water Resources Research 57 (2): e2020WR027463. https://doi.org/10.1029/2020WR027463

Somers, L.D.; McKenzie, J.M. 2020. A review of groundwater in high mountain environments. WIREs Water 7 (6): e1475. https://doi.org/10.1002/wat2.1475

Uhlemann, S.; Dafflon, B.; Peterson, J.; Ulrich, C.; Shirley, I.; Michail, S.; Hubbard, S.S. 2021. Geophysical monitoring shows that spatial heterogeneity in thermohydrological dynamics reshapes a transitional permafrost system. Geophysical Research Letters 48 (6): e2020GL091149. https://doi.org/10.1029/2020GL091149

Villarroel, C.D.; Ortiz, D.A.; Forte, A.P.; Tamburini Beliveau, G.; Ponce, D.; Imhof, A.; López, A. 2022. Internal structure of a large, complex rock glacier and its significance in hydrological and dynamic behavior: A case study in the semi-arid Andes of Argentina. Permafrost and Periglacial Processes 33: 78-95. https://doi.org/10.1002/ppp.2132

Vittecoq, B.; Reninger, P.-A.; Lacquement, F.; Martelet, G.; Violette, S. 2019. Hydrogeological conceptual model of andesitic watersheds revealed by high-resolution heliborne geophysics. Hydrology and Earth System Sciences 23 (5): 2321-2338. https://doi.org/10.5194/hess-23-2321-2019

Walvoord, M.A.; Kurylyk, B.L. 2016. Hydrologic impacts of thawing permafrost-a review. Vadose Zone Journal 15 (6): 1-20. https://doi.org/10.2136/vzj2016.01.0010

Wang, J.-Z.; Jiang, X.-W.; Wan, L.; Wang, X.-S.; Li, H. 2014. An analytical study on groundwater flow in drainage basins with horizontal wells. Hydrogeology Journal 22: 1625-1638. https://doi.org/10.1007/s10040-014-1146-9

Zhu, L.; Wu, S.; Zhang, C.; Misra, S.; Zhou, X.; Cai, J. 2023. Characterization of pore electrical conductivity in porous media by weakly conductive and nonconductive pores. Surveys in Geophysics 44: 877-923. https://doi.org/10.1007/s10712-022-09761-w

Zinn, B.A.; Konikow, L.F. 2007. Potential effects of regional pumpage on groundwater age distribution. Water Resources Research 43 (6): W06418. https://doi.org/10.1029/2006WR004865