Spatiotemporal variation of modern lake, stream, and soil water isotopes in Iceland

As global warming progresses, changes in high-latitude precipitation are expected to impart long-lasting impacts on Earth systems, including glacier mass balance and ecosystem structures. Reconstructing past changes in high-latitude precipitation and hydroclimate from networks of continuous lake records offers one way to improve forecasts of precipitation and precipitation–evaporation balances, but these reconstructions are currently hindered by the incomplete understanding of controls on lake and soil water isotopes. Here, we study the distribution of modern water isotopes in Icelandic lakes, streams, and surface soils collected in 2002, 2003, 2004, 2014, 2019, and 2020 to understand the geographic, geomorphic, and environmental controls on their regional and interannual variability. We find that lake water isotopes in open-basin (through-flowing) lakes reflect local precipitation, with biases toward the cold season, particularly in lakes with sub-annual residence times. Closed-basin lakes have water isotope and deuterium excess values consistent with evaporative enrichment. Interannual and seasonal variabilities of lake water isotopes at repeatedly sampled sites are consistent with instrumental records of winter snowfall; summer relative humidity; and atmospheric circulation patterns, such as the North Atlantic Oscillation. Summer surface soil water isotopes span the entire range of seasonal precipitation values in Iceland and appear to be consistently overprinted by evaporative enrichment, which can occur throughout the year, although the sampling depths were shallower than rooting depths for many plant types. This dataset provides new insight into the functionality of water isotopes in Icelandic environments and offers renewed possibilities for optimized site selection and proxy interpretation in future paleohydrological studies on this North Atlantic outpost.

1 Introduction

The Arctic hydrological cycle is directly tied to broad Earth systems, such as the cryosphere and biosphere, and is expected to change rapidly with continued anthropogenic warming (Serreze et al., 2006). By the end of the century, Arctic precipitation in some regions is forecasted to increase by up to 40 % – with a progressively larger proportion of rain – due to sea ice retreat, increased local evaporation, and enhanced poleward moisture transport (Bintanja and Selten, 2014; Bintanja and Andry, 2017; Bintanja et al., 2020). As a result, long-lasting impacts on hydrology, glacier mass balances, ocean circulation, ecosystem structures, and global climate system feedbacks are expected (Douville et al., 2021). However, regional-scale predictions of future precipitation changes feature large uncertainties spatially and across models owing to large decadal-scale natural variability and inter-model differences (Bintanja et al., 2020). Developing records that quantify changes in Arctic precipitation beyond the instrumental period in Earth's history provides one means of quantifying how precipitation has changed with past climate variability and can lead to more refined predictions of future precipitation and its impact on other Earth systems (e.g., Linderholm et al., 2018; Konecky et al., 2023). Geologic records of past precipitation are, however, similarly limited at present because there are few ubiquitous proxies that clearly relate to specific aspects of hydrologic cycling (Sundqvist et al., 2014).

Lake sediments are a commonly used archive to reconstruct past environmental change due to their wide geographical distribution, sensitivity to local and regional perturbations, and continuous integration of environmental information at the catchment scale (Adrian et al., 2009). In the Arctic, various physical (e.g., grain size, sediment accumulation rates, and varves), biological (e.g., pollen, chironomids, and diatoms), and geochemical proxies (e.g., compound-specific stable isotopes) have been used to infer past patterns of precipitation from lake sediment (Linderholm et al., 2018). Of these, the hydrogen or oxygen stable-isotopic compositions of aquatic and terrestrial macrofossils and lipids are commonly used because they are abundant in lakes, and a dominant proportion of the H or O in the preserved structure is derived from plant growth water, which can then be used to estimate past water isotopes through time (Verbruggen et al., 2011; Rach et al., 2017; Holtvoeth et al., 2019). Changes in the past isotopic composition of precipitation, as recorded by water isotope proxies, can therefore yield important information on the response of regional hydrologic cycles to global changes in climate.

However, the stable-isotope signatures preserved in sedimentary archives must be carefully considered with respect to the growth habitat and reservoir of meteoric water that source organisms are likely to use, changes in the seasonal bias of precipitation, catchment geometry, and changes in landscape dynamics (Sachse et al., 2012; Cluett and Thomas, 2020; McFarlin et al., 2023). Efforts to map the spatiotemporal variability of modern lake water isotopes in highly seasonal climates have advanced our understanding of the processes that control the isotopic composition of lake water on short and long timescales (Gibson et al., 2002; Leng and Anderson, 2003; Jonsson et al., 2009; Cluett and Thomas, 2020; Gorbey et al., 2022; Kjellman et al., 2022), yet such mapping efforts have not yet been conducted in Iceland. Furthermore, while regional studies are not uncommon for lake water isotopes, similar systematic efforts are limited for soil water isotopes, meaning there is a gap in observations for the growth water utilized by terrestrial plants as it isotopically relates to precipitation. Most calibration work for compound-specific water isotope proxies (e.g., n-alkyl lipids) uses mean annual precipitation (e.g., Sachse et al., 2012), yet there are local observations that demonstrate that Arctic soil water can be seasonally biased (Cooper et al., 1991; Throckmorton et al., 2016; Bush et al., 2017; McFarlin et al., 2019; Muhic et al., 2023) and/or evaporatively enriched (Tetzlaff et al., 2018; Eensalu et al., 2023). It is also common to assume that terrestrial plant growth water is restricted to seasonal precipitation because Arctic plant growth is limited to late spring and summer (e.g., Balascio et al., 2013; Curtin et al., 2019; Kjellman et al., 2020; Thomas et al., 2020). These environmental distinctions are important because lake water (aquatic) and precipitation (terrestrial) proxies are often used comparatively to reconstruct hydrologic budgets through time, and differences in the isotopic composition of varying compounds have been used quantitatively to estimate changes in relative humidity (Rach et al., 2017) or seasonal precipitation amounts (Corcoran et al., 2021; Thomas et al, 2020).

In this study, we present modern lake, stream, and soil water δ²H and δ¹⁸O data from a network of sites across Iceland (Fig. 1a). We compare lake and soil water isotopes collected during the summer and winter field seasons in 2002, 2003, 2004, 2014, 2019, and 2020, with climatological, geomorphic, and geographic parameters to assess which factors may account for the observed isotopic variability. We also present a nearly year-long time series of lake outflow water isotopes collected in 2002. Our data show that, in general, lake water isotopes (δ²H_lake) in Iceland reflect cold-season or annual precipitation with an overprint of evaporative enrichment in closed-lake basins. The interannual variability of δ²H_lake is then related to the specific hydrologic budget (i.e., winters with differing amounts of snowfall and/or years of low relative humidity) and atmospheric circulation patterns for individual years. In contrast to lake water, surface soil water isotopes (δ²H_soil) broadly demonstrate values that align with summer-biased precipitation and annual evaporative enrichment despite the fact that most precipitation is delivered in the autumn and winter. Collectively, these data have important implications for the interpretation of water isotope proxies in Icelandic lake sediment records and provide an objective framework within which to estimate past precipitation isotopes in local Holocene archives.

Figure 1 Overview map of Iceland and modern climate. (a) Key lake and soil locations from 2014, 2019, and 2020 are shown with blue circles and are numerated (see Table 1 for site information), and additional lakes, ponds, and streams from 2002, 2003, and 2004 are shown in black (see the Supplement for site information). (b–e) Mean and standard deviation of monthly climate and water isotope data from Reykjavík, Iceland, for the years 1992–2018: (b) mean air temperature (°C) and (c) mean precipitation (mm) from the IMO (IMO, 2023) and (d) δ²H (‰, VSMOW) and (e) d-excess values (‰, VSMOW) from GNIP (IAEA/WMO, 2015).

Table 1 Key sample site information.

NA denotes “not available”, i.e., the lake water depth was not taken to calculate lake volume and residence time.

2 Regional setting

Iceland (∼ 103 000 km²) is situated in the middle of the northern North Atlantic Ocean. Due to its location at the confluence of major oceanic and atmospheric patterns, Iceland's terrestrial climate is sensitive to variations in the strength and position of these systems (e.g., Hanna et al., 2006; Geirsdóttir et al., 2013, 2020). The southern and western coastlines of Iceland are bathed in relatively warm, saline Atlantic waters, whereas the northern and, occasionally, the eastern coastlines are influenced by cool and relatively fresh Arctic and polar waters (Stefánsson, 1962). Iceland is also proximal to the Icelandic Low, a semi-permanent center of low atmospheric pressure that forms one dipole of the North Atlantic Oscillation (NAO), along with Azores High, and is associated with frequent cyclone activity. This pressure difference between the north and south is the dominant mechanism that controls the strength and location of westerly wind and storm tracks across the northern North Atlantic (Hurrell et al., 2003). Given the different oceanographic regimes surrounding Iceland, incoming air masses are strongly influenced by the Atlantic or Arctic waters over which they originate.

The modern climate of Iceland is generally characterized by strong winds, frequent precipitation, mild winters, and cool summers (Einarsson, 1984; Hanna et al., 2004). The coldest part of Iceland is the central highlands, and snowfall in winter is more common in the northern regions. The south coast of the island is warmer, wetter, and windier compared to the northern coast (Hanna et al., 2004, 2006). In Reykjavík, southwest Iceland (Fig. 1a), temperatures average near freezing in the winter months and 10 °C in summer months (Fig. 1b, IMO, 2023). On monthly timescales, the precipitation amount in Reykjavík is inversely related to temperature, with a greater proportion of precipitation falling during autumn and winter compared to during summer (Fig. 1c, IMO, 2023). However, on longer-term decadal timescales, the precipitation amount shows a moderate positive correlation with temperature (Hanna et al., 2004). The isotopic composition of precipitation recorded at the Global Network of Isotopes in Precipitation (GNIP) station in Reykjavík is, on average, most depleted during the winter months and most enriched in the month of May, although there is large year-to-year variability (Fig. 1d, IAEA/WMO, 2015). This variability reflects the competing moisture sources to Iceland and likely explains a lack of correlation between δ²H of precipitation and air temperature (Fig. S1a in the Supplement). Deuterium excess (d-excess) values of precipitation are generally highest during the winter and lowest during the summer (Fig. 1e, IAEA/WMO, 2015). Reykjavík`s annual local meteoric water line (LMWL, Fig. 2, IAEA/WMO, 2015) of

$$\begin{array}{}\text{(1)}& {\mathit{\delta }}^{\mathrm{2}}\mathrm{H}=\mathrm{7.78}\times {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}+\mathrm{5.29}\end{array}$$

has a similar slope but lower d-excess values than the global meteoric water line (GMWL, Craig, 1961) of

$$\begin{array}{}\text{(2)}& {\mathit{\delta }}^{\mathrm{2}}\mathrm{H}=\mathrm{8}\times {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}+\mathrm{10}.\end{array}$$

LMWLs of the different seasons (winter, DJF; spring, MAM; summer, JJA; and autumn, SON) are variable, where summer plots distinctly below the annual LMWL, while winter, spring, and autumn plot closer to and typically above the LMWL (Fig. S1b). In hydrologically closed-lake surface water systems, kinetic fractionation from evaporation will result in isotopically enriched water relative to precipitation and low d-excess values (Feng et al., 2016; Bowen et al., 2018). These are more appropriately described by a local evaporation line (LEL, Fig. 2), which plots below the LMWL (Gat, 1996).

Figure 2 Lake water stable isotopes from Iceland plotted against the GMWL (Craig, 1961), LMWL (GNIP Reykjavík, IAEA/WMO, 2015), and LEL (closed lakes, this study). Samples are separated by (a) hydrologic connectivity for all lakes (open and closed, n= 223); (b) year and season of collection for all lakes (n= 223); (c) residence time for open lakes from 2014, 2019, and 2020 (RT, n= 61); and (d) d-excess values for all lakes (n= 223). Panels (e) and (f) show annual LMWLs (open lakes, E) and LELs and slopes (closed lakes, f), respectively, based on lake water isotope values in comparison to the LMWL derived from GNIP data in Reykjavík (black line). Raw GNIP data are plotted as gray points.

The mountainous topography of the island, particularly in the southeast, produces a strong orographic effect on precipitation, with annual precipitation in excess of 5000 mm atop glacial summits, such as Vatnajökull and Mýrdalsjökull (Einarsson, 1984; Crochet et al., 2007). Along with southerly to southeasterly prevailing winds, precipitation amount tends to be heaviest in the south and southeast of Iceland, with a rain shadow over the central highlands, especially north of Vatnajökull (Einarsson, 1984). This geography results in distillation of air masses as they move inland and, thus, progressively more depleted δ²H values of precipitation at higher altitudes farther from the coastline (Bödvarsson, 1962; Friedman et al., 1963; Árnason, 1976). In addition, temperature and variations in moisture source also exert strong controls on precipitation isotopes in Iceland (Bödvarsson, 1962; Friedman et al., 1963) and across the North Atlantic (Dansgaard, 1964; Jouzel et al., 2000; Bowen and Revenaugh, 2003; Sodemann, et al., 2008; Steen-Larsen et al., 2015). Moisture sources in Iceland are particularly susceptible to variations of the NAO, which impacts the strength and location of westerly winds and storm tracks on interannual to decadal cycles, particularly during winter months (Hurrell et al., 2003; Olsen et al., 2012). The impact of the NAO is apparent in the modern δ²H and d-excess values of precipitation in Reykjavík, which show negative ($r=-\mathrm{0.37}$) and positive correlations (r=0.45) with the NAO index, respectively (Baldini et al., 2008).

3 Methods and materials

3.1 Sample collection

We present a compilation of water isotope data collected over multiple field campaigns in 2002, 2003, 2004, 2014, 2019, and 2020 (Fig. 1a). Between 2002 and 2004, we collected lake water samples from 135 lakes and ponds (mainly in summers, with several sites being revisited in February and March; see Fig. 2b), as well as geothermal streams (n= 2), non-geothermal streams and lake water outflow (n= 8), and snow meltwater (n= 3) (Fig. S3). The sample set from 2003 (n= 148) is particularly important as it contains a large number of closed-lake basins that allow us to make robust estimates of the LEL (Fig. 2a–b). The sample sets from 2002 (n= 20) and 2004 (n= 3) are comparatively small (Fig. 2b), and since the same lakes were not revisited between the 3 years, the 2002–2004 datasets are not used to assess the detailed interannual variability of lake water isotopes. For one lake in NW Iceland (Efstadalsvatn), Ragna Aðalsteinsdóttir, the farmer residing at the site, collected water samples from the lake's outflow every 2 weeks between March and November 2002. We then visited a total of 26 lakes in summer 2014 (n= 14), summer 2019 (n= 17), winter 2020 (n= 4), and summer 2020 (n= 19) (Fig. 2b and Table 1). Of the lakes visited during the summer field season in 2014, seven were revisited in summer 2019 and 2020. All lakes visited during the summer field season in 2019 were also visited during the subsequent summer field season in 2020. This repeat-sampling strategy for 2014–2020 lake water forms the basis of our analysis focused on the interannual variability of lake water isotopes. During each field season, we sampled surface water from the center of each lake using Nalgene bottles. For select lakes in summer 2019 (n= 2) and winter 2020 (n= 4), we also sampled bottom water from the center of the lake using a Niskin water sampler. All water samples were immediately transferred to glass vials with no headspace, sealed with Teflon and/or electrical tape to prevent evaporation, and refrigerated prior to isotopic analysis. Further details on sampling locations, additional proxy data, and associated water chemistry for the 2014, 2019, and 2020 datasets can be found in Florian (2016), J. H. Raberg et al. (2021) and J. Raberg et al.(2021b), and Raberg et al. (2023).

We sampled a total of 23 surface soil samples from 18 locations in summer 2019 (Fig. 1a and Table 1). Of the 18 locations, 16 are within the catchments of lakes sampled for lake water during the same field season as described above. The two remaining soil sites (Skriðdalur and Mývatn) do not have corresponding lake water samples as they were targeted to provide forested endmembers. Soil samples were collected to evenly integrate the top 10 cm of soil. Each soil sample was sampled in 50 mL Falcon tubes, which were fully filled, sealed with Teflon tape, and kept refrigerated to prevent evaporation of water before isotopic analysis. Further details on sampling locations and additional proxy data datasets can be found in Raberg et al. (2024).

3.2 Water isotope analysis

Lake water samples were analyzed in the Stable Isotope Laboratory at the Institute of Arctic and Alpine Research, University of Colorado Boulder. Samples from 2002 to 2004 were measured via mass spectrometry methods using the CO₂ equilibration method for oxygen isotopes (Epstein and Mayeda, 1953) and the uranium reduction method for hydrogen isotopes (Vaughn et al., 1998). For the remainder of the sample set, stable-isotope measurements (δ²H and δ¹⁸O) were conducted using a Picarro L2130-i cavity ring-down mass spectrometer. Each sample was run six times, with the first three data points being discarded to prevent memory isotope effects. Mean values were calculated using the final three measurements. A replicate was run after ∼20 samples to check for reproducibility. We report data using δ notation in per mil (‰) relative to the Vienna Standard Mean Ocean Water (VSMOW):

$$\begin{array}{}\text{(3)}& \mathit{\delta }=\mathrm{1000}\times \left({\displaystyle \frac{R_{\mathrm{sample}}}{R_{\mathrm{VSMOW}}}}-\mathrm{1}\right),\end{array}$$

where R is ²H $/$ ¹H and ¹⁸O $/$ ¹⁶O. Samples were analyzed against a suite of secondary laboratory standards that have been vigorously calibrated to primary reference materials, namely VSMOW2 and SLAP. The secondary laboratory standards (Florida water, Boulder water, and Antarctic water) have δ²H $/$ δ¹⁸O values of $-\mathrm{2.82}/-\mathrm{0.76}$ ‰, $-\mathrm{111.65}/-\mathrm{14.15}$ ‰, and $-\mathrm{239.13}/-\mathrm{30.30}$ ‰, and an analytical error of 1 ‰ and 0.1 ‰, respectively.

All soil samples were sent to the University of Utah Stable Isotope Research Facility for soil water vacuum line extraction and isotope measurements (see Cook et al. (2017) for methodology). Values are reported in delta notation relative to VSMOW.

Deuterium excess (d-excess) values were calculated after Dansgaard (1964) as follows:

$$\begin{array}{}\text{(4)}& \text{d-excess}={\mathit{\delta }}^{\mathrm{2}}\mathrm{H}-\mathrm{8}\times {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}.\end{array}$$

Given that δ²H and δ¹⁸O values are generally correlative, we focus our discussion on δ²H and d-excess values for simplicity.

3.3 Additional site information

For each location, we compiled and computed additional information to test which environmental parameters are related to lake and soil water δ²H values.

For all lakes, we determined whether they were open or closed basins via the presence or absence of surface outflow in summer 2023 satellite imagery, and we measured elevation (m) and the closest distance from coast (km) using Google Earth (2023). We recognize there may be some ambiguity in this binary open–closed designation as some lakes that are, at times, open basins can be closed in late summer or in a year of low snowpack, for example. For lakes visited in 2014, 2019, and 2020 only, water depths were measured using a weighted measuring tape from the lake depocenter, and lake surface area (km²) and lake catchment area (km²) were measured in QGIS using the ÍslandsDEM (2×2 m resolution, Landmælingar Íslands). A complete model of lake bathymetries is unavailable, and so we estimate lake volumes (V, km³) for the 2014–2020 lakes based on half of an ellipsoid:

$$\begin{array}{}\text{(5)}& V={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\left({\displaystyle \frac{\mathrm{4}}{\mathrm{3}}}\times A_{\mathrm{lake}}\times Z\right),\end{array}$$

where A_lake is the lake surface area, and Z is the maximum depth of the lake. Lake water residence times (RTs) were calculated after Jonsson et al. (2009):

$$\begin{array}{}\text{(6)}& \mathrm{RT}={\displaystyle \frac{\mathrm{lake}\mathrm{volume}}{\mathrm{runoff}}},\end{array}$$

where

$$\begin{array}{}\text{(7)}& \mathrm{runoff}=\mathrm{catchment}\mathrm{area}\times \mathrm{precipitation}\times {\displaystyle \frac{R}{P}},\end{array}$$

with $R/P$ being the runoff $/$ precipitation ratio (Gibson and Edwards, 2002), where the latter is less than 1 because some precipitation that falls in a lake catchment is lost to evaporation, transpiration, sublimation, or groundwater flow. Quantitative information on these processes for Icelandic lakes is limited, and so we use an $R/P$ value of 0.5, which is consistent with existing Arctic lake studies conducted elsewhere (Gibson et al., 2002; Gorbey et al., 2022). We acknowledge that the $R/P$ value may vary spatially between sites and between years for individual lakes experiencing changes in the precipitation–evaporation balance. Additionally, the role of groundwater is an important unknown. However, as we are more broadly interested in the correlation of residence time with water stable isotopes, an estimate of 0.5 is deemed to be sufficient to address this variable. Finally, we obtained mean annual air temperature (MAT), mean annual precipitation (MAP), and seasonal means from the weather station closest to lakes sampled in 2014–2020, operated by the Icelandic Meteorological Office (IMO, 2023), including the number of years preceding sample collection equivalent to individual lake residence times. As the elevation of the closest weather station often differs from the elevation of the lake site, we also scaled temperature integrations for each lake according to the global lapse rate (6 °C km⁻¹).

For soil sites, we used Google-Earth-derived information on the elevation and distance from the coastline, as described above for lake sites. In addition, we measured soil water content (SWC),

$$\begin{array}{}\text{(8)}& \mathrm{SWC}={\displaystyle \frac{\mathrm{wet}\mathrm{mass}-\mathrm{dry}\mathrm{mass}}{\mathrm{dry}\mathrm{mass}}}\times \mathrm{100},\end{array}$$

and in situ soil temperatures for one year using iButton loggers deployed at 10 cm soil depth in summer 2019 and collected in summer 2020 (J. Raberg et al., 2021). iButton temperatures were recorded every 3 h, which we then used to produce mean annual and seasonally integrated values. We additionally used iButtons to measure in situ surface and bottom water temperatures at 6 h intervals over the same period for four lake sites (J. Raberg at al., 2021). Finally, we include information on the dominant vegetation cover above each soil sample that we observed in the field during collection: forest (Betula pubescens); dwarf shrub (Betula nana); heath shrub (Ericaceae); moss, lichen, and/or prostrate; grass; and moss. Except for Torfavatn, which is classified as a histosol (TOC > 20 %), all collected soils are classified as brown-gleyic andosols (TOC < 12 %) (Arnalds, 2004).

3.4 Statistical analyses

We computed Pearson correlation matrices and p values to assess the relationship between water isotopes, d-excess values, and various environmental parameters using the open-source software R (R Core Team, 2022) and the corrplot package (Wei and Simko, 2021). To assess the significance of lake and soil water isotope population differences and of water isotopes and environmental variable correlations, we performed Student's t tests. We define significant results to be those with p values < 0.05.

4 Results

For all lakes (2002–2020), the δ²H of lake water (δ²H_lake) ranges from −89.53 ‰ to −1.35 ‰ and generally falls within the natural range of meteoric water isotopes recorded at the Reykjavík GNIP station (Fig. 2). Spatially, δ²H_lake values tend to be more depleted in northern Iceland compared to in the south (Fig. S2). The d-excess values of lake water range from −11 ‰ to 15.71 ‰, where the values of open lakes are typically higher, ranging from −4.30 ‰ to 15.71 ‰ (Fig. 2d), and exhibit no clear spatial pattern (Fig. S2). The δ²H values of stream, snowmelt, and geothermal spring water values also fall within the natural range of meteoric water isotopes, ranging from −89.55 ‰ to −55.66 ‰, −79.11 ‰ to −60.13 ‰, and −79.37 ‰ to −73.93 ‰, respectively (Fig. S3), and the d-excess values range from −6.58 ‰ to 14.91 ‰. In general, lake water isotopes for all open lakes plot along the LMWL, whereas water isotopes for closed lakes evolve along a LEL (Fig. 2a). Using data from all closed lakes, the LEL is described by the following relationship:

$$\begin{array}{}\text{(9)}& {\mathit{\delta }}^{\mathrm{2}}\mathrm{H}=\mathrm{5.90}\times {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}-\mathrm{13.09}(r=\mathrm{0.97},p\mathit{<}\mathrm{0.05}).\end{array}$$

Using data from 2014, 2019, and 2020, open lakes that have residence times of less than 1 year (efficient flushing) follow the annual LMWL closely. In contrast, open lakes that have residence times of over 1 year (less efficient flushing) feature larger deviations from the LMWL (Fig. 2c).

Figure 3 Pearson correlation matrices for (a) all open and closed lakes and (b) only open-lake water δ²H (‰, VSMOW); d-excess values; and various climatological, geographic, and environmental parameters. Shown are elevation and distance from the coastline (n= 223); lake depth, catchment area, and lake volume (n= 80); residence time and annual and seasonal precipitation amount (n= 61); and annual and seasonal temperature (n= 80). MAP refers to mean annual precipitation, MWP refers to mean winter precipitation, MFP refers to mean autumn precipitation, MSpP refers to mean spring precipitation, MSP refers to mean summer precipitation, MAT refers to mean annual temperature, MWT refers to mean winter temperature, MFT refers to mean autumn temperature, MSpT refers to mean spring temperature, MST refers to mean summer temperature. Empty boxes reflect statistically insignificant correlations (p > 0.05). Note that lake waters sampled in 2002–2004 are only used in the regression analyses of the elevation and distance from the coast.

Figure 4 Interannual variability of lake surface water δ²H compared with climatological datasets. (a) δ²H (‰, VSMOW) of surface water from seven lakes collected in July 2014, February 2019, and September 2020; (b) corresponding d-excess values (‰, VSMOW); (c) mean monthly air temperature (°C) for each lake; (d) mean monthly precipitation (mm) for four lakes (not available for Leirvogsvatn, Margrétarvatn, and Högnavatn); (e) mean monthly relative humidity (%) for four lakes (not available for Heiðarvatn, Margrétarvatn, and Högnavatn); and (f) mean monthly NAO (NOAA National Centers for Environmental Information, 2023). For (f), NAO annual means are denoted by dashed horizontal lines and are calculated based on the 12 months prior to lake water sampling (red indicates positive, and blue indicates negative).

For all 2014–2020 lakes (open and closed), we find strong correlations between δ²H_lake and mean annual (r= 0.66, p < 0.05), mean winter (r= 0.67, p < 0.05), mean autumn (r= 0.69, p < 0.05), and mean spring temperatures (r= 0.64, p < 0.05) and a weak but still significant correlation with mean summer temperature (r= 0.38, p < 0.05) (Fig. 3a). However, we find no significant correlations between mean annual precipitation amount and δ²H_lake values, nor do we find any correlation with seasonal precipitation amount (Fig. 3a), which may be in part due to the inherent spatial variability of precipitation and the use of weather stations with varying distances to lakes. If closed-lake basins are excluded, we find similar relationships between δ²H_lake and temperature (Fig. 3b). However, excluding closed lakes also reveals that open lakes are weakly correlated with mean spring (r= 0.37, p < 0.05) and mean summer precipitation amounts (r= 0.33, p < 0.05) (Fig. 3b). Using the entire 2002–2020 lake dataset, we find a moderate correlation between δ²H_lake values and elevation (r= -0.46, p < 0.05) and a weak correlation between δ²H_lake values and the distance from the coastline ($r=-\mathrm{0.33}$, p < 0.05) (Fig. 3). Lake water d-excess values show weak to moderate correlations with most variables in all lakes, including elevation, the distance from the coast, precipitation amount, and temperature (Fig. 3a). However, when closed lakes are excluded, d-excess values are no longer correlated with precipitation amounts, except for mean autumn precipitation ($r=-\mathrm{0.44}$, p < 0.05) (Fig. 3b). For repeatedly sampled lakes, summer δ²H_lake values are relatively more enriched in 2019 than summer δ²H_lake values in 2014 and 2020 (Figs. 4a and 5). Summer d-excess values are also relatively lower in 2019 compared to in 2014 and 2020 (Figs. 4b and 5). For the four lakes that have seasonal data from 2019 and 2020, winter surface δ²H_lake values are more depleted than summer surface δ²H_lake values, with differences ranging from −0.54 % to 16.51 % (Fig. 6a). Winter lake water d-excess values are relatively higher than summer d-excess values, except for one lake, Stóra Viðarvatn (Fig. 6b). For the 2002 lake outflow stream time series from Efstadalsvatn (Fig. 1a), δ²H values range from −78.75 ‰ at the end of April to −70.68 ‰ at the end of June and are generally ²H-enriched between June and August (Fig. 7a), and d-excess values range from 11.13 ‰ to 14.29 ‰, with the highest values being obtained in June (Fig. 7b).

Figure 5 Difference in September 2020 and 2019 lake surface water δ²H (‰, VSMOW) vs. d-excess values (‰, VSMOW) from 17 lakes. Squares denote lakes located on the NW highlands (Fig. 1), and black and gray fillings denote lakes with residence times > 1 and < 1 year, respectively.

Figure 6 Interannual variability of lake surface and bottom water δ²H compared with climatological datasets. (a) δ²H (‰, VSMOW) of surface and bottom water from four lakes collected in September 2019, February 2020, and September 2020; (b) corresponding d-excess values (‰, VSMOW); (c) difference between in situ lake surface and bottom water temperature (°C) (not available for Litla Viðarvatn; J. Raberg et al., 2021); (d) mean monthly air temperature (°C) for each lake; and (e) mean monthly precipitation (mm) for each lake (note that Litla Viðarvatn is the same as Stóra Viðarvatn).

Figure 7 Time series of (a) stream outflow δ²H (‰, VSMOW) and (b) d-excess values (‰, VSMOW) from Efstadalsvatn between March and October 2002 in comparison to local instrumental weather data: (c) air temperature (°C) and (d) precipitation (mm). Reykjavík GNIP precipitation δ²H (‰, VSMOW) is shown in panel (a) (gray) for comparison (IAEA/WMO, 2015).

The δ²H values of surface soil water (δ²H_soil) from 2019 range from −91.34 ‰ to −31.71 ‰ and fall within the natural range of meteoric water isotopes recorded at the Reykjavík GNIP station across all seasons (Fig. 8). Due to the small size of the δ²H_soil dataset compared to lakes, it is difficult to assess if there is any clear spatial variability (Fig. S4). Soil water d-excess values range from −5.04 ‰ to 9.84 ‰ (Fig. 8b). In general, surface soil water isotopes fall below the LMWL, although samples collected from beneath densely vegetated moss or grass cover plot along the LMWL (Fig. 8). Five of the sites have two soil samples from the same lake catchment, with the difference in δ²H_soil values ranging from 1.11 ‰ to 21.65 ‰. In general, surface soil waters are ²H-enriched compared to surface lake waters from the same year, although four soil water samples are ²H-depleted compared to surface lake water at the same locations. For temperature, we find a moderate correlation between δ²H_soil and in situ soil mean summer temperature (r= 0.51, p < 0.05) and no significant correlations between mean annual, mean winter, mean autumn, and mean spring temperature and δ²H_soil values (Fig. 9). For precipitation, we find a moderate correlation between δ²H_soil and mean summer precipitation amount (r= 0.56, p < 0.05) and no significant correlation between mean annual, mean winter, mean autumn, and mean spring precipitation amount and δ²H_soil values (Fig. 9). In contrast, soil water d-excess values are positively correlated with mean annual (r= 0.58, p < 0.05), mean winter (r= 0.61, p < 0.05), mean autumn (r= 0.56, p < 0.05), and mean spring precipitation amounts (r= 0.50, p < 0.05) but not with the mean summer precipitation amount (Fig. 9). We also find weak to moderate correlations between d-excess values and mean annual (r= 0.58, p < 0.05), mean autumn (r= 0.57, p < 0.05), mean spring (r= 0.50, p < 0.05), and mean summer precipitation amounts (r= 0.33, p < 0.05) (Fig. 9). We find moderate correlations between δ²H_soil and elevation ($r=-\mathrm{0.43}$, p < 0.05) and between d-excess values and elevation (r= 0.43, p < 0.05) and no significant correlation with the distance from the coast (Fig. 9). Finally, we find no significant correlation between δ²H_soil and soil water content (Fig. 9).

Figure 8 Soil water stable isotopes (n= 23) plotted against the GMWL and LMWL. (a) Samples colored according to the overlying vegetation as described during collection in September 2019 and (b) samples colored by d-excess values (‰, VSMOW).

Figure 9 Pearson correlation matrix of 2019 soil water δ²H (‰, VSMOW); d-excess values; and various climatological, geographic, and environmental parameters (n= 23). Empty boxes reflect statistically insignificant correlations (p > 0.05).

Student's t tests for lake and soil water populations show that δ²H and δ¹⁸O values are not statistically different between the two groups (p > 0.05). However, Student's t tests for d-excess values, which take in account both δ²H and δ¹⁸O variables, show that lake and soil water populations are statistically different from one another (p < 0.05). This supports the separate assessment of lake and surface soil water isotopes in Iceland.

5 Discussion

5.1 Lake water isotopes

5.1.1 Environmental controls

Icelandic lake water isotopes from 2002 to 2020 generally plot close to or on the LMWL (estimated using GNIP data from Reykjavík, IAEA/WMO, 2015), demonstrating that most lake water in Iceland is isotopically consistent with unmodified precipitation and a low-temperature, high-precipitation environment. Lakes with δ²H_lake values that generally plot on the LMWL are all open lakes (n= 149, Fig. 2a) with active inflows providing recharge from meteoric water that falls in the catchment and is stored according to lake residence times (Jonsson et al., 2009). For 2014 to 2020 lakes, the strongest predictor of δ²H_lake in open lakes is cold-season temperatures: winter (r= 0.53, p < 0.05), autumn (r= 0.62, p < 0.05), and spring (r= 0.63, p < 0.05) (Fig. 3b). However, due to the lack of correlation between GNIP δ²H values of precipitation and air temperature (Fig. S1a), we do not interpret δ²H_lake as a function of temperature but rather a function of precipitation sourced from cool seasons. This is supported by the fact that open lakes with residence times of less than 1 year (efficient flushing) plot close to the LMWL (Fig. 2c), which was found to be biased towards the cold season in Reykjavík GNIP precipitation data (Fig. S1b), and the fact that a larger proportion of precipitation is delivered in autumn and winter (Fig. 1c). Surprisingly, we observe no relationship between precipitation amount and δ²H_lake (Fig. 3b), which may be due to the large standard deviation of seasonal precipitation amounts recorded at the Reykjavík GNIP station (Fig. 1c). Open lakes that have residence times of greater than 1 year (less efficient flushing) deviate from the LMWL more than lakes with residence times of less than 1 year (Figs. 2c and S1). Following the seasonal bias of LMWLs (Fig. S1b), this indicates that the dominant proportion of precipitation falling in these lake catchments falls in the cold season (plots close to annual LMWL) and possibly in summer (plots below annual LMWL). For the latter “summer” lakes, the low corresponding d-excess values imply that there is some degree of seasonal evaporation (Fig. 2d), meaning that precipitation may still be dominantly sourced from cold seasons but that repeated episodes of evaporation may obscure the seasonal signal in lakes with longer residence times.

Several other mechanisms also help to explain the trends observed in Icelandic lake water isotopes. First, 49 lakes in the entire dataset (2002 to 2020) are closed basins, meaning that the meteoric water that falls in the catchment is contained until it is lost through evaporation. In Iceland, the effect of evaporation (kinetic fractionation) is observed through the deviation of closed lakes' water isotope values along the LEL and relatively low d-excess values (Fig. 2d). Given the maritime and humid climate of Iceland, as well as the fact that some of the evaporative lakes are located on the coast, this observation contrasts with documented lake water isotopes in the higher-latitude Arctic, where evaporation is minimal in coastal lakes (Leng and Anderson, 2003; Cluett and Thomas, 2020; Kjellman et al., 2022). Compared to these other Arctic lake studies, we hypothesize that the relatively warmer temperature, low amounts of summer precipitation, and high wind speeds in Iceland, particularly along the coast and in the highlands (Einarsson, 1984), may contribute to enhanced seasonal evaporation during the relatively longer open-water season in these lakes (e.g., Feng et al., 2016; Gonfiantini et al., 2020). In addition, as winds are typically strongest in winter (Einarsson, 1984), a partially ice-free cool season, when the water–air temperature difference is greatest, may contribute to enhanced lake water evaporation as well. Second, we observe moderate and weak correlations in all lakes (2002 to 2020) between δ²H_lake and two geographic parameters, namely elevation (r= −0.46, p < 0.05) and the distance from the coastline (r= −0.33, p < 0.05) (Fig. 3a). Icelandic lakes that are further inland tend to be at higher elevations; thus, Rayleigh distillation is also a control on the spatial variability of lake water isotopes in Iceland, consistently with other surface water, near-surface groundwater, and precipitation stable-isotope studies (Bödvarsson, 1962; Friedman et al., 1963; Árnason, 1976; Stefánsson et al., 2019).

5.1.2 Interannual variability

LMWLs calculated from open lakes in the years 2002, 2003, 2014, 2019, and 2020 reveal similarities and differences that are likely to be driven by broad annual changes in climate across Iceland (Fig. 2e). First, for open lakes, we observe that 2002 stands out with a LMWL slightly above the other years, which are all similar in slope and intercept (Fig. 2e). The samples from 2002 were all collected from the northwestern peninsula, which is characterized by a cooler climate and a greater proportion of precipitation falling as snow (Einarsson, 1984). Following our seasonal LMWL analysis of Reykjavík GNIP data, the higher intercept of the 2002 lake water LMWL relative to the GNIP LMWL is consistent with a greater proportion of isotopically depleted winter precipitation (Fig. S1b). The lower intercepts of lake water LMWLs across the other years (2003, 2014, 2019, and 2020; Fig. 2e) are consistent with spring and autumn LMWLs (Fig. S1b) and likely reflect the inclusion of lakes that receive relatively lower amounts of winter snow. In any case, these data demonstrate that, broadly, the seasonality of lake water isotopes is consistent with cold-season precipitation.

LELs calculated from closed lakes in the years 2003, 2014, 2019, and 2020 also reveal similarities and differences that are likely to be driven by broad changes in Icelandic climate (Fig. 2f). The intersection of the LEL, defined based on closed-lake basin water isotopes, and the LMWL of precipitation is commonly used to infer the isotopic composition of lake water input prior to evaporative enrichment (e.g., Gibson et al., 2016). The annual LELs show that 2019 has the highest LEL–LMWL intercept, followed by 2020, 2014, and 2020, all of which fall within the range of GNIP data from Reykjavík (Fig. 2f). Some high-latitude studies infer such changes to reflect changes in the seasonality of precipitation, where higher (lower) LEL–LMWL intercepts reflect more summer (winter) precipitation (e.g., Kjellman et al., 2022). However, the seasonality of Reykjavík precipitation spans the range of LEL–LMWL intercepts and, as previously noted, is better reflected in seasonal LMWLs (Fig. S1b). However, we do note that the slope of the annual LELs derived from lake water isotopes varies between years, where 2019 is the shallowest (m= 5.5), followed by 2020 (m= 5.8), 2003 (m= 5.9), and 2014 (m= 6.7) (Fig. 2f). Under the assumption that LELs are predominately the result of variations in the precipitation–evaporation balance and that all lakes respond similarly, the slope of LELs can be used as an indicator of relative humidity variability in high-latitude lakes, where lower slopes reflect increased lake water evaporation (e.g., Leng and Anderson, 2003; Kopec et al., 2018). Hence, these data suggest that 2019 was relatively dry and that 2014 relatively wet compared to 2003 and 2020.

Seven lakes in Iceland were sampled for surface lake water isotopes during September 2014, 2019, and 2020 and demonstrate interannual δ²H_lake variability ranging from 0.24 (Margrétarvatn, NW Iceland) to 17.88 ‰ (Heiðarvatn, E Iceland) (Fig. 4a) and d-excess variability ranging from 1.69 (Torfdalsvatn, N Iceland) to 10.83 (Högnavatn, NW Iceland). In all seven lakes, δ²H_lake values in 2019 were relatively enriched and had lower d-excess values than in 2014 and 2020. For some lakes, winter snowfall may be a factor explaining δ²H_lake values, particularly for those at higher elevations. At Heiðarvatn (442 m a.s.l., E Iceland), ²H-depleted snow from the preceding winter may have yielded a winter-biased signal that lingered into the summer season for some years. We note that the winter preceding September 2014 featured substantially higher amounts of winter snowfall (1063 mm) than the winter preceding September 2019 (553 mm) at Heiðarvatn (Fig. 4d), which is consistent with the relatively ²H-depleted values in 2014 (17.88 ‰, Fig. 4a). For the low-elevation, coastal lake Torfdalsvatn (N Iceland), the amount of winter precipitation did not vary substantially between 2014 and 2020, as it did for Heiðarvatn. However, local relative humidity at Torfdalsvatn was significantly lower through the winter, spring, and summer in 2019 compared to in 2014 and 2020 (Fig. 4e). As Torfdalsvatn is a closed lake, the relative enrichment of δ²H_lake values and the lower d-excess values in 2019 (Fig. 4a–b) are consistent with the relatively dryer conditions experienced during the same year.

A total of 17 lakes were sampled in September of both 2019 and 2020 and revealed interannual variability, with δ²H_lake value differences between 2019 and 2020 ranging from −3.31 (Skjaldarvatn, W Iceland) to 15.89 ‰ (Margrétarvatn, NW Iceland) and d-excess value differences ranging from -1.70 (Skjaldarvatn, W Iceland) to 5.65 ‰ (Berufjarðarvatn, NW Iceland) (Fig. 5). In general, we observe lower δ²H_lake values and higher d-excess values in 2020 compared to in 2019, except for one lake (Skjaldarvatn), which has higher δ²H_lake values and lower d-excess values in 2020 (Fig. 5). We predicted that lakes with residence times of less than 1 year (efficient flushing) would feature larger interannual variability (e.g., Gibson et al., 2002; Leng and Anderson, 2003). However, lake residence time does not correlate with the magnitude of the δ²H_lake value difference between 2019 and 2020 (Fig. 5). However, we do observe greater δ²H_lake and d-excess difference for lakes in the NW highlands (squares) compared to for lakes elsewhere in Iceland (circles, Fig. 5). We note that, during the 2020 sample collection, there was fresh snow accumulation in the catchments of the NW highland lakes in September, whereas we did not observe this for the same time period in 2019 (Fig. S5). This is consistent with more ²H-depleted values in 2020 relative to 2019; thus, snowmelt at high elevations may partially explain the relatively ²H-depleted values for NW highland lakes compared to lower-elevation sites elsewhere in Iceland.

Finally, given that lake water isotope values generally reflect unmodified precipitation (Fig. 2), the precipitation source also likely impacts the Icelandic δ²H_lake values. Variation of the NAO will impact the strength and location of westerly winds and storm tracks over interannual to decadal cycles. In NAO+ (strong low pressure over Iceland), the westerlies center over Iceland, resulting in warm and moist air being delivered from the south. Conversely, during NAO− (weak low pressure over Iceland), the storm track shifts southward, resulting in relatively cooler and drier conditions in Iceland (Hurrell et al., 2003). Therefore, we would expect that, during NAO+ (NAO−) years, lake water in Iceland may be ²H-enriched (²H-depleted) in through-flowing lakes due to different source regions from the south and north, respectively, which can differ isotopically by 100 ‰ (Bowen and Revenaugh, 2003; Hurrell et al., 2003). In closed-lake basins, we may expect ²H enrichment and lower d-excess values of lake water in NAO− years due to the drier conditions that promote evaporation. In 2019, the average NAO index for the 12 months preceding lake water sample collection was negative (−0.21), whereas those for 2014 and 2020 were positive (0.23 and 0.09, respectively; Fig. 4f, NOAA National Centers for Environmental Information, 2023). Our dataset demonstrates more ²H-enriched water isotopes with lower d-excess values in 2019 relative to 2014 and 2020 (Fig. 4a–b). For through-flowing lakes, the ²H enrichment in 2019 may be the result of less ²H-depleted winter snowfall in a NAO− year, whereas, for closed lakes (e.g., Torfdalsvatn), the ²H enrichment and low d-excess values may additionally reflect evaporative enrichment resulting from the changes in relative humidity induced by NAO− conditions. For the latter, this is also consistent with the relatively lower LEL slope in 2019 (Fig. 2f), which has previously been shown to correlate with NAO variability in lakes in western Greenland (Kopec et al., 2018). In any case, annual changes in the NAO index and the associated changes in moisture source and relative humidity are consistent with the isotopic variability of lake water in Iceland, demonstrating that the precipitation source is likely a contributing influence in relation to δ²H_lake values across Iceland on at least interannual timescales.

5.1.3 Seasonal variability

Four lakes were sampled for surface lake water three times in two different seasons between 2019 and 2020: September 2019, February 2020, and September 2020. In addition, we took bottom water samples from two lakes in September 2019 (Torfdalsvatn and Heiðarvatn) and from all four lakes in February 2020 (Torfdalsvatn, 4.4 m depth; Heiðarvatn, 13 m depth; Stóra Viðarvatn, 19.8 m depth; and Litla Viðarvatn, 1.5 m depth) (Fig. 6a). In general, we observe more depleted surface δ²H_lake values in winter, which likely reflects the seasonal change in precipitation, where winter snow is isotopically depleted relative to rain, which dominates the other seasons (Fig. 6a). In addition, winter d-excess values of surface water are higher than in summer, except in Stóra Viðarvatn, which may reflect varying degrees of summer evaporation in closed (Torfdalsvatn) and even open lakes (Heiðarvatn and Litla Viðarvatn). Compared to surface lake water, winter bottom δ²H_lake values are relatively depleted in Torfdalsvatn, relatively enriched in Heiðarvatn and Litla Viðarvatn, and equivalent in Stóra Viðarvatn (Fig. 6a). Winter bottom water d-excess values are relative higher in Torfdalsvatn and Litla Viðarvatn, lower in Heiðarvatn, and similar in Stóra Viðarvatn (Fig. 6b). While the difference varies between lakes, the general trends are likely to reflect the winter stratification of the water columns. Using in situ daily lake water temperatures from the same depths that δ²H_lake and d-excess values are from for three of the lakes (Torfdalsvatn, Heiðarvatn, and Stóra Viðarvatn; Fig. 6c, J. Raberg et al., 2021), we observe strong winter water column stratification in Torfdalsvatn and Heiðarvatn but not in Stóra Viðarvatn, which is consistent with differences in surface and bottom δ²H_lake and d-excess values, where Torfdalsvatn and Heiðarvatn's δ²H_lake values differ by 4.92 ‰ and 8.62 ‰, respectively, while Stóra Viðarvatn's values differ only by 0.06 ‰ (Fig. 6a). For Litla Viðarvatn, in the absence of continuous bottom water temperature measurements, sonde water quality measurements (i.e., lake water temperature, pH, dissolved oxygen, and conductivity) from September 2019 and February 2020 (Raberg et al., 2023) demonstrate that this lake is also stratified during the winter (Fig. S6), which is consistent with the difference between surface and bottom water δ²H_lake values of 2.74 ‰ (Fig. 6a).

In addition to seasonal δ²H_lake values, we also collected bimonthly water samples from the outflow of a low-elevation lake in northwestern Iceland (Efstadalsvatn, Fig. 7a, black) between March and November 2002. Efstadalsvatn has a short residence time of ∼0.1 per year (Table 1), meaning that meteoric water should be efficiently recharged through the lake. Compared to the seasonal changes in precipitation δ²H values in Reykjavík for the same year (Fig. 1d), the δ²H values of the outflow water are relatively dampened, ranging from −78.75 ‰ at the end of April to −70.68 ‰ at the end of June (Fig. 7a, black), but these are similar in scale compared to the interannual and seasonal variability observed in δ²H_lake values from other lakes (Figs. 4–6). Notably, we observe relatively ²H-depleted values between March and April, followed by rapid (14 d) ²H enrichment of 5.83 ‰ at the beginning of June and then rapid (2 d) ²H depletion of 7.56 ‰ at the end of August, after which δ²H values are stable and low (Fig. 7a, black). The relatively ²H-depleted values in March and April may be due to the seasonal variability of precipitation δ²H values in Reykjavík (Fig. 7a, gray) and/or lingering snowmelt in the catchment associated with relatively cool seasonal air temperatures (Fig. 7c). The ²H enrichment observed at the beginning of June occurs at the same time as an increase in the local precipitation amount (Fig. 7d), which is likely to be in the form of rain (as opposed to ²H-depleted snow) due to above-freezing air temperatures (Fig. 7c). Evaporative enrichment is not likely to explain this step change as Efstadalsvatn has an open catchment and a short residence time, and d-excess values are high and inconsistent with a kinetic fractionation effect (Fig. 7b). This ²H enrichment of lake outflow water is also superimposed over a more gradual ²H enrichment of Reykjavík precipitation values (Fig. 7a). The cause of the rapid ²H depletion at the end of August may be related to the source and/or seasonal variability of precipitation δ²H values. In Reykjavík, monthly averages show a large depletion between August and September (Fig. 7a, gray). As precipitation sources may differ between Reykjavík (southwest) and Efstadalsvatn (northwest), these similarities are difficult to confidently correlate. In any case, while there are no data for Efstadalsvatn, surface and bottom water isotopes and temperatures from Torfdalsvatn, a similarly shallow lake (3.8 m vs. 5.8 m depth, respectively; Caseldine et al., 2003), show that it is dimictic and seasonally stratified from early April through October (Fig. 6a–c). This suggests that Efstadalsvatn could also be seasonally stratified, resulting in August rainfall flushing through the lake surface without mixing with the deeper water (e.g., Boehrer and Schultze, 2008).

5.2 Surface soil water isotopes

Isotopes of surface soil water collected in summer in Iceland span the entire range of seasonal precipitation values for Iceland, suggesting that the seasonality of source water and/or the residence time of soil water are spatially variable across the island. Values correlate significantly with mean summer environmental parameters, i.e., precipitation amount (r= 0.56, p < 0.05) and temperature (r= 0.52, p < 0.05) (Fig. 9), suggesting that summer is the source water season for at least some sites and that there is, broadly, more integration of warm-season precipitation than of lake water across the soil water dataset. In contrast, d-excess values correlate with mean annual (r= 0.58, p < 0.05), mean winter (r= 0.61, p < 0.05), mean autumn (r= 0.56, p < 0.05), and mean spring precipitation amounts (r= 0.50, p < 0.05) (Fig. 9), demonstrating that the input of seasonal precipitation to soils in Iceland is not clear-cut. As most Icelandic surface soil water isotopes plot below the LMWL and have low to negative d-excess values (Fig. 8), the pattern of soil waters isotopes, regardless of source water seasonality, is broadly consistent with modification via evaporation. While the residence time of Icelandic soil water is unknown, the residence time of soil water in other wet, low-energy North Atlantic locations is typically controlled by soil type and overlying vegetation (Tetzlaff et al., 2014; Geris et al., 2015a, b; Sprenger et al., 2017). The andic properties of the soil types in Iceland (i.e., histosol and brown-gleyic Andosol) result in relatively high water retention (Arnalds, 2004), which is consistent with the correlations of precipitation and d-excess values throughout the year. As soil temperatures at these sites periodically dip below freezing during winter months (Fig. S7), seasonal freezing of the ground may contribute to the integration of summer precipitation δ²H_soil values, despite lower summer precipitation amounts, by preventing the incorporation of cold-season water into a frozen-water–soil matrix at some sites. Alternatively, surface soil water residence times at some sites could just be short enough that, by summer sampling, winter precipitation has percolated deeper into the soil profile. When soil samples are separated by the type of vegetation cover, we observe that moss- and grass-covered soil samples associated with marsh environments plot closer to the LMWL and have higher d-excess values compared to other types of vegetation cover (Fig. 8). This likely reflects the saturated marsh environment that these soils were sampled from, which have high water tables, retain water for longer than soils with efficient drainage, and are less prone to evaporation.

For five locations, we took two surface soil samples from the same lake catchments but different locations, as measured for lake water isotopes (Heiðarvatn, Másvatn, Stóra Viðarvatn, Systravatn, and Torfdalsvatn). The δ²H_soil value differences for the samples from the same lake catchments range from 1.11 (Systravatn) to 21.65 ‰ (Heiðarvatn), and corresponding d-excess value differences range from 0.54 (Systravatn) to 5.66 (Heiðarvatn). The variability observed in δ²H_soil and d-excess values from the same location suggests that additional processes beyond climate and geography contribute to soil water isotopic variability in Iceland. This heterogeneity may be one reason why correlations are weaker between δ²H_soil values and mean summer temperature (r= 0.52, p < 0.05), precipitation amount (r= 0.52, p < 0.05), and elevation ($r=-\mathrm{0.43}$) (Fig. 9) compared to the same correlations found with δ²H_lake values (Fig. 3). While vegetation cover may explain some spatial variability in the broader Icelandic dataset, it likely does not account for intersite variability. For example, the sites with high intersample δ²H_soil value differences (Heiðarvatn, with 21.65 ‰, and Stóra Viðarvatn, with, 11.37 ‰) are both covered by the same vegetation type (moss and heath shrub, respectively), while the only location that had different vegetation cover (Másvatn, grass and heath shrub) had one of the lowest δ²H_soil value differences (2.04 ‰). Variations in snow cover duration may be another possible explanation for inter-site variability. However, soil temperatures (3 h resolution) for the following year show that sites with larger δ²H_soil value differences (e.g., Heiðarvatn and Stóra Viðarvatn) have melt-out differences of 3 and 8 d, respectively, whereas sites with small δ²H_soil value differences (e.g., Másvatn, Systravatn, and Torfdalsvatn) show melt-out differences with a larger range from 3 to 13 d (Fig. S7). Therefore, the parameters for which we have data, i.e., local climate parameters, fine-scale observations of vegetation cover, and in situ thermal regimes, do not adequately explain the site-to-site variability of surface soil water isotopes.

5.3 Implications for paleoclimate proxies

Lake water isotope proxies have become an increasingly common tool to infer past changes in the high-latitude hydrologic cycle, including the seasonality of precipitation (e.g., Kjellman et al., 2020; Thomas et al., 2020; Corcoran et al., 2021), moisture balance (e.g., Anderson and Leng, 2004; Balascio et al., 2018), and precipitation source (e.g., Cowling et al., 2022; van der Bilt et al., 2022). Foundational to any paleohydrological study is a clear understanding of the local environmental controls behind water isotope variability. For Icelandic lake water, we find that cold-season precipitation, lingering snowpack in fringe seasons, summer relative humidity, and NAO cyclicity are important modulators of δ²H_lake on interannual timescales, particularly for lakes with sub-annual residence times. These are important considerations as they indicate that the seasonality of runoff may not always align with the seasonality of precipitation. Combined with the differences we observe between surface and bottom water δ²H_lake, the complexities of δ²H_lake are consistent with the poorer relationships observed between typically aquatic mid-chain plant waxes and δ²H_lake values vs. typically terrestrial long-chain plant waxes and δ²H of precipitation values (McFarlin et al., 2019). Surface soil water isotopes are highly variable across the island, with only some sites being supportive of previous assumptions that terrestrial biomarker water isotopes reflect summer precipitation (e.g., Balascio et al., 2013; Curtin et al., 2019; Kjellman et al., 2020; Thomas et al., 2020). However, we stress that these soil water isotopes are from the surface (i.e., upper 10 cm of the soil), and evidence suggests that high-latitude plant roots may source water from the surface (Amin et al., 2020), as well as from below the surface layer (e.g., Eensalu et al., 2023). Therefore, these data may or may not reflect broader plant growth water, and the enrichment observed in Icelandic surface soils also may not reflect water found in deeper soils. Hence, while still valuable, the soil data currently provide a limited perspective on the soil water isotope systematics needed to inform plant growth and terrestrial leaf wax water isotope studies in Iceland.

In Iceland, there are only two Holocene records of paleo-precipitation developed from water isotopes proxies, one using aquatic invertebrates (i.e., chironomids) from a lake in northeastern Iceland (Stóra Viðarvatn, Wooller et al., 2008) and another using leaf wax n-alkanes from a fjord in northwestern Iceland (Moossen et al., 2015). Both records observe inconsistencies with known temperature variability and, through multi-proxy comparisons, infer that the changes in water isotopes were related to precipitation. The dominant controls that we find in modern lake water isotopes support prior assumptions that chironomid δ¹⁸O values from Stóra Viðarvatn encode information about past changes in precipitation seasonality and moisture source (Wooller et al., 2008). However, we do note that Stóra Viðarvatn has one of the longest residence times in the dataset (∼17 years, Table 1), and the lake's modern water isotopes appear to be insensitive to interannual climate variability (Fig. 6). While residence time is not a strong predictor of lake water isotope variability in Iceland, the relative insensitivity of Stóra Viðarvatn to annual changes in moisture supply likely results in a decadally averaged paleohydrological record. On the other hand, the fjord record from northwestern Iceland is based on terrestrial leaf waxes that were subsequently transported to the fjord (Moossen et al., 2015), which are complicated by unknown transport and residence times, as well as by uncertainties in marine age model reservoir corrections. The relatively limited perspective that our modern surface soil water isotopes provide regarding terrestrial biomarker water isotope proxies also makes it challenging to interpret the fjord record, although it remains likely that the fjord leaf wax δ²H values are summer biased. Looking forward, the strong linear relationship between lake water and aquatic-invertebrate (e.g., chironomids) oxygen isotopes (Verbruggen et al., 2011) suggests that their application in Icelandic lake sediment records is promising. However, further calibration of terrestrial and aquatic plant waxes, their isotopic variability, their relationship to source waters, and their relative contribution to the sedimentary pool is vital to further advance the use of biomarker-derived water isotope proxies in Icelandic sedimentary records (e.g., McFarlin et al., 2019; Dion-Kirschner et al., 2020; Hollister et al., 2022).

Our modern lake and soil water isotope dataset provides important insights into the impact that various geographic, geomorphic, and environmental variables have on the distribution of water isotopes across Iceland. Importantly, the data open up the possibility for optimized site selection and more informed paleohydrological studies in Iceland. For example, we identify lakes that are sensitive to evaporative enrichment, whose sediments can be useful in reconstructing Holocene moisture balances, and lakes that record unmodified meteoric water, which offers the possibility of reconstructing changes in precipitation seasonality and moisture source related to atmospheric patterns. Collectively, these paleohydrological records can help address open questions related to, for example, the impact of precipitation on the Holocene mass balance of Icelandic ice caps and shifts in terrestrial plant communities. As both glaciers and terrestrial plants are part of positive feedback loops in the climate system through changes in albedo and evapotranspiration (Serreze and Barry, 2011), understanding the paleohydrology of Iceland is critical to better understand the past evolution of climate, glaciers, and ecology and to predict future climate and environmental scenarios over the 21st century.

6 Conclusions and recommendations

We report Icelandic lake, stream, and soil water isotopes (δ²H and δ¹⁸O) measured from samples collected in 2002, 2003, 2004, 2014, 2019, and 2020. Guided by statistical correlation analyses between δ²H; d-excess values; and a suite of climatological, geomorphic, and geographic parameters, we draw the following conclusions.

Lake water isotopes reflect local precipitation, with biases toward the cold season, particularly for lakes with sub-annual residence times, which is consistent with the time of year when most precipitation is delivered. A subset of lakes, which are closed basins, demonstrate δ²H_lake and d-excess values affected by evaporative enrichment. This evaporation likely occurs during the typically open-water summer season when relative humidity is lower, lake water surface temperatures are elevated, and the lake surface is exposed to wind. Lake water isotopes are generally ²H-depleted at higher elevations and at locations farther from the coast, consistently with Rayleigh distillation.

Interannual variabilities of lake water isotopes and d-excess values are consistent with observations that 2019 received less winter snowfall, was relatively dry during the summer, and was characterized by a negative NAO, all contributing to a lower LEL slope and ²H-enriched and low-d-excess lake water relative to 2014 and 2020.

Seasonal variability of lake water isotopes shows that lake surface water is relatively ²H-depleted, with higher d-excess values in February 2019 compared to in September 2019 and 2020, likely due to the predominance of snow over rain in winter months. Additionally, bottom water isotopes reflect winter stratification that is independently recorded in the seasonal temperature and water chemistry cycles of the lakes.

Isotopes and d-excess values of summer surface soil water reflect local precipitation overprinted by evaporative enrichment that can occur throughout the year. The soil water values span the seasonal range of precipitation, suggesting heterogeneity in soil water seasonality. Water isotope measurements from samples collected just tens of meters apart can show large variability, suggesting that additional factors beyond the local climate, especially topography and vegetation, which are inherently related, may be important factors.

The collective dataset optimizes site selection for future studies in Icelandic paleohydrology focusing on changes in moisture balance, precipitation seasonality, and moisture source and provides a process-based framework to guide downcore water isotope proxy interpretations. Specifically, sedimentary records from closed-basin lakes offer the potential to reconstruct Holocene changes in the precipitation–evaporation balance, while sedimentary records from open-basin lakes, particularly those with sub-annual residence times, may provide insight into past changes in moisture source related to winter snowfall variability and atmospheric circulation patterns. At present, the simplest sedimentary records of precipitation would be derived from aquatic-invertebrate (e.g., chironomids) oxygen isotopes given the strong linear relationship observed between lake water and aquatic-invertebrate oxygen isotopes. Future work is needed to constrain the relationship between aquatic plant waxes and their isotopic variability and their relationship to lake water before secure paleohydrological studies can be developed from them. Similarly, given the complexity and uncertainty observed in Icelandic surface soil water isotopes, we urge caution in the application of terrestrial water isotope proxies, such as long-chain n-alkyl lipids, before further calibration studies are conducted.