Late Summer Water Sources in Rivers and Lakes of the Upper Yana River Basin, Northern Eurasia, Inferred from Hydrological Tracer Data

Abstract

Major ions, stable isotopes, and trace elements, including rare earth elements (REEs), are used as natural tracers in the qualitative assessment of potential water sources in lakes and rivers of the upper Yana River basin, between Verkhoyansk and Chersky Ranges, during the late summer period. Three distinct regions were sampled, and a dominant water source in each region was qualitatively inferred from water chemistry data. The REE distribution pattern was found to be highly regional and controlled by pH and carbonate contents. Mountain headwater stream at the Verkhoyansk Range north slope, the Dulgalakh River, shows an input from a mixture of shallow groundwater and icing meltwater, with a depleted isotopic signature (δ¹⁸O below –21‰), d-excess (d_ex = δ2H − 8·δ18O) above 18, enrichment in Mg and Sr, and depletion in heavy REEs. The Derbeke Depression lakes and streams are fed by rainfall having ultra-low total dissolved solids (TDS) content, below 25 mg/L, and a convex-up REE pattern. In a medium mountainous river at the Chersky Range flank, the Dogdo River, leaching through fissured Jurassic carbonates is a dominant runoff pathway. Riverine water is heavily depleted in light REEs, but enriched in Mo, Rb, Sb, W and U. In the Dulgalakh River water, high positive Sm and Gd anomalies were observed, attributed either to local geology (greenshists), historical mining legacy, or contemporary winter road operations.

1. Introduction

Frozen soils of permafrost regions are an important stock of chemical elements currently withdrawn from geochemical cycling [1], and contain over 800 Pg C of soil organic carbon [2,3]. Ground ice and ice-rich Ice Complex terrain also store about 43 TgC as organic carbon (DOC) and 33.6 TgC as inorganic carbon (DIC), in dissolved phase vulnerable to thawing [4]. Arctic climate warming and projected permafrost degradation are expected to release substantial amounts of organic carbon and chemical elements into the hydrographic network, reintroducing them into global turnover [5,6,7,8]. Shorter water transfer pathways (i.e., melting massive ground ice in exposures along river banks), or physical permafrost disturbances assure immediate delivery of dissolved material to rivers and streams [9].

Permafrost environments host a wide variety of hydrological processes affecting the regime of both surficial water bodies and groundwater aquifers [10,11,12,13]. Soil moisture redistribution and seasonal migration of the active layer base, transient water storage in solid phase, and baseflow redistribution across timescales through such storage are important features of permafrost hydrology, driving solute transport as well [11,14]. In continuous permafrost, water transfer is confined to the active layer and local open (through) taliks under lakes and major rivers [13,15,16]. In discontinuous permafrost, hydrological connectivity extends to multiple non-frozen zones within catchments, allowing better surface-subsurface connectivity and increased groundwater drainage to streams [17,18,19,20]. Wildfires [21,22] and higher snow accumulation, owing to climate change or snow retention in shrubs and snow fences [23,24,25], lead to thermally-driven permafrost degradation, and development of non-merging permafrost, also termed ‘isolated talik’ [26,27], ‘lateral talik’ [28], ‘residual talik’, or ‘residual thaw layer’ [18,29,30], ‘perennial thaw zone’ [31]. Non-merging permafrost conveys water toward streams and rivers during winter, when other streamflow sources are inactive, leading to increased winter daily flows across the global Arctic [32,33,34].

Hydrological tracers serve to differentiate between varying water sources, highlighting runoff pathways, and groundwater discharge and residence time in permafrost-dominated catchments [35,36]. In various settings, basic water chemistry, stable water isotopes, dissolved organic carbon, iron and rare earth elements, and organic compounds were monitored to study water transfer through complex permafrost terrain in a single- or multi-proxy analysis. Stable water isotopes are widely used to infer water origin, travel times and dominant pathways in surface runoff and evapotranspiration across permafrost regions [37,38,39,40,41]. In an extensive field study along the permafrost gradient in western Siberia, dissolved organic carbon, major, and trace elements reacted to latitudinal change in permafrost continuity, and were found useful in estimating the seasonal variability in input from peatlands to the hydrological system [42,43,44]. At the sub-seasonal scale, basic water chemistry (electrical conductivity, pH) and DOC observations during a high-frequency sampling campaign revealed seasonality in DOC origin and pathways in a medium catchment of the Northern Yenisey region [45]. In permafrost groundwater studies, radioactive tritium, radon (²²²Rn), and chlorofluorocarbons, which only recently appeared in the atmosphere, are used jointly to infer groundwater recharge and residence times in intra- and supra-permafrost aquifers [46,47,48,49]. Tritium release from thawing permafrost may significantly enrich surficial waters and groundwater seepage, especially in discontinuous permafrost [50].

Permafrost thaw is expected to significantly alter water chemistry [51,52,53]. Thus, in the long-term monitoring studies, the latter may serve a reliable proxy for permafrost degradation processes. Thawing permafrost directly releases meltwater to the subsurface compartment, changing pre-thaw water chemistry [54,55]. In the Canadian High Arctic, an abrupt increase in sulphate ion content was attributed to thermally-driven permafrost degradation [56]. At the same time, after forest fires and permafrost thaw, DOC concentrations were found to decrease with deeper active layer and longer pathways in mineral soils [44,57].

Increased groundwater discharge caused by thawing permafrost is another cause of changes in water chemistry [28,58,59,60]. Groundwater icings, or ‘aufeis’, at the interface between deep subsurface and surface compartments, withdraw water from winter base flow and redistribute it toward summer and autumn [61]. Their presence affects hydrological regime and fluvial patterns in mountainous permafrost regions of Russia [62,63], northern Canada [64,65], and Alaska [66]. With ongoing climate warming, icings disappear significantly earlier in summer, reflecting not only higher air temperatures, but also less icing volume, and fewer icings are expected to develop in the warmer future [67,68]. This implies changes to stream chemistry as more deep groundwater is expected to enter streams during winter, but less icing meltwater will be available, and thus allows tracking the effects of permafrost degradation. Differentiating between the baseflow and icing meltwater by natural tracer methods may be a challenge, as both originate from the same groundwater reservoirs. Hence, specific natural tracers accounting for the transformation of chemical composition in freeze-defreeze cycles are needed to ‘fingerprint’ icing discharge in rivers.

This study presents the first data on basic water chemistry, stable water isotopes (δ²H and δ¹⁸O), and trace elements, including rare earth elements (REEs), in streams and lakes of a remote and isolated inland region in the Verkhoyansk and Chersky Range system, Northern Eurasia. The hydrology of mountainous headwater streams is controlled by permafrost through groundwater routing and is expected to respond rapidly to changing climate and permafrost thawing [69]. The main objectives of this study were the exploratory description of water chemistry in this remote locality, for the first time in the literature, and qualitative assessment of dominating water sources in the sampled water bodies using natural tracer data.

2. Study Area

Helicopter survey was carried out in the north-eastern part of Sakha (Yakutia) Republic, Northern Eurasia (Figure 1). The region is remote and extremely hard to access, and sparsely populated, with population density below 0.1/km². Regional climate is harsh continental: it is one of the coldest locations in the Northern Hemisphere, with the lowest observed air temperature below −70 °C and mean annual air temperature as low as −14 °C (

Table 1. Mean monthly and annual air temperature, Verkhoyansk Meteo station, 1991–2020.

Jan.	Feb.	Mar.	Apr.	May	Jun.	Jul.	Aug.	Sep.	Oct.	Nov.	Dec.	MAAT
−44.7	−42.2	−28.7	−11.0	4.3	13.8	16.4	12.2	2.6	−13.5	−33.9	−43.4	−14.0

). Mean annual precipitation is between 250 and 300 mm. In the last few decades, the study region experiences pronounced warming that particularly affected summer temperatures, showing statistically significant upward trends [70]. In 2020, Verkhoyansk Meteo station recorded +38.0 °C, recognized by World Meteorological Organization in 2021 as a new Arctic temperature record [71].

The surveyed region is located in the upper section of the Yana River basin; sampling campaign covered three distinct locations in the headwaters of its major tributaries: (1) the northern slope of the Verkhoyansk Range, in the headwaters of Dulgalakh R. and Nelgese R.; (2) the north-western flank of the Elga Highlands, known as Derbeke Depression, at the interfluve of the Adycha and Derbeke Rivers; and (3) the south-western slope of the Chersky Range, the Tuostakh River basin (Figure 1). Region 1 is a piedmont region with altitudes from 1300 to 1800 m a.s.l., dissected by deeply incised glacier valleys with N-S orientation. Region 2 is presumably a Pleistocene glacier terminus, an elevated plain with altitudes from 900 to 1000 m a.s.l., and abundant lakes. Region 3 is a mountain chain with altitudes from 1800 to 2300 m a.s.l.

Sampled rivers are typical mountainous rivers, with wide anabranching channels and numerous icing glades. They mostly have pluvial regime with low spring freshet owing to minor snow accumulation during winter, and multiple major rain events during short summer. Hydrological change in the region is observed in autumn and winter months, and represents a step shift rather than a monotonic trend, while summer flows remain relatively stable. Increased May discharge may be attributed to higher icing meltwater runoff [72].

Regional geology and geomorphology were previously studied in the scope of geological prospections [73], as the region is rich in gold and associated Sb and As, Sn, W, Co, Pb, Zn, Bi, In, Cd, rare earth elements (REEs). In Regions 1 and 2, pre-Quaternary basal rocks are Triassic argillites and conglomerates, in Region 3, Neogene sands of the lower Dogdo River reach change upstream to Jurassic sandstones, and heavily folded and cut sandstones and dolomites with ages from Ordovician to Devonian.

Continuous permafrost in the region in 250 to 500 m deep, with well-developed periglacial relief, including tabular ground ice, pingos, ice wedges and byllars, alas depressions and thermokarst lakes in river valleys and on elevated flatlands, hillslope water tracks on gentle rolling slopes, and kurums and blockfields at water divides [74]. Major river valleys follow fault lines and tectonic depressions, allowing intense subpermafrost groundwater discharge. Subsequently, icings and former icing glades are abundant in most river valleys of the region. Some icing remnants were observed during the sampling campaign on the Dulgalakh and Nelgese Rivers, Verkhoyansk Range, and their tributaries (see Figure 1), and they are also known to be present also in the Tuostakh River basin, including the sampled Dogdo River. Minor glaciers were studied in the Yana River basin, but their total volume and area are negligible [72], and limit their contribution to river runoff in the basin. Former glaciations have left their imprint in local topography, by creating circular moraine features in piedmont areas.

3. Materials and Methods

The helicopter-based sampling campaign took place from 22 to 24 August 2017. The helicopter was in lease for non-scientific purposes, and thus the flight plan was prepared and approved regardless of research activities. Owing to limited resources, sampling strategy was opportunistic and only covered water bodies close to helicopter landing sites selected by the helicopter crew. In rivers and lakes, samples were taken at a depth of 10 to 20 cm from the water surface, with hands protected by latex powder-free Nitrile gloves. Samples for major ions (n = 7) were collected in 0.25 L high-density polyethylene (HDPE bottles), samples for trace elements (n = 7) and stable water isotopes (n = 12), collected in 15 mL HDPE sterile conical centrifuge tubes. All samples were sealed with Parafilm tape and stored in the field at temperatures between 5 °C and 10 °C in cool place, in organic topsoil, and in the laboratory, at 4 °C, before the analysis. Samples for trace elements were pre-filtered through 0.22 μm cellulose acetate membrane filters before analysis.

Major ion content, electric conductivity, and pH, were measured at the in-house facilities of Permafrost Groundwater and Geochemistry Laboratory, P.I. Melnikov Permafrost Institute, SB RAS (Yakutsk, Russia), following Russian national guidelines. Hydrocarbonates content was recalculated from total alkalinity, measured by titration with 20% accuracy [75]. Major anions and cations were measured by capillary ion electrophoresis [76,77], with accuracy ca. 15%, using ‘Kapel’ capillary ion electrophoresis system (Lumex, Russia). Analytical results were expressed in mg/L and milliequivalents per liter (meq/L), then anion and cation equivalent masses were summed up and equaled 100% to calculate %-equivalent masses (%-eq.) for each anion and cation.

Stable water isotopes were measured using multiflow isotope-ratio mass spectrometry (IRMS) at SHIVA platform, Laboratory of Functional Ecology, UMR5245 CNRS-UT3-INP (Toulouse, France). The analytical precision of the method is ±0.1‰ for δ¹⁸O, and ±1.0‰ for δ²H. Each water sample was measured in duplicate and averaged, so each δ¹⁸O/δ²H value presented in the paper is a mean value. Deuterium excess (d-excess, or d_ex, ‰) was calculated as d_ex = δ²H − 8·δ¹⁸O. Trace elements, including REEs, were measured at the Institute of Microelectronics Technology and High Purity Materials, RAS (Chernogolovka, Russia) using ICP-MS/AES, at ng/L precision. REEs concentrations were NASC-normalized using values from [78,79].

Statistical analysis was performed in RStudio [80], a graphical interface for R language [81]. Hellinger-transformed absolute element contents were used to calculate Bray-Curtis dissimilarity using vegdist() function from the package ‘vegan‘ [82]. Agglomerative clustering of water samples was done using agnes() function from the package ‘cluster’ [83], and elements clustering was done using fviz_cluster() function from the package ‘factoextra’ [84].

4. Results

4.1. Major Ions

Sampled waters mostly have ultra-low total dissolved solids (TDS) content, below 100 mg/L. Lowest TDS values were observed in creeks and lakes of Region 2 (Derbeke Depression), and ranged from 14.7 to 21.3 mg/L. The Dulgalakh River (Region 1) and the Dogdo River (Region 3) exhibited elevated TDS contents of 50.7 mg/L, and 105 mg/L, respectively. For anions (in meq/L), hydrocarbonate ion dominates at 91 to 99% in Region 2 samples, but is slightly less abundant, 73–74%, in the Dulgalakh and Dogdo River waters, where the sulphate-ion is also present (Figure 2). The dominant cation (in meq/L) is Ca²⁺ at 49 to 69% followed by Mg²⁺ at 30 to 45% (Figure 2 and

Table 3. Major ions content in samples from the Yana River headwaters (descriptions and location, see Table 2); numerator shows concentrations in mg/L, denominator, in %-eq. The numerator is above the underline and the denominator is below it.

Lab #	pH	Eh, mV	Ca2+	Mg2+	Na+	K+	NH4+	HCO3−	SO42−	Cl−	NO2−	NO3−	Hard 1	TDS
373/17B	6.71	410	12.1	4.09	0.6	0.04	0.1	42.5	11.5	0.42	0.01	0.2	0.942	50.8
			61.6	34.2	2.6	1.04	0.56	73.2	25.2	1.24	0.02	0.34
374/17B	7.71	342	27.0	7.37	0.4	0.8	0.2	88.4	24.7	0.42	0.01	0.112	1.952	105.2
			67.3	30.3	0.84	1.02	0.54	73.3	26.0	0.60	0.01	0.09
375/17B	6.12	493	2.7	1.88	0.05	0.3	0.1	18.9	0.02	0.04	0.01	0.112	0.289	14.7
			44.2	50.8	0.70	2.50	1.80	98.9	0.10	0.37	0.06	0.57
376/17B	6.07	459	3.51	1.97	0.05	0.2	0.2	21.3	0.016	0.104	0.01	0.2	0.337	16.9
			49.3	45.5	0.62	1.44	3.14	98.1	0.10	0.83	0.06	0.91
377/17B	6.16	515	5.94	1.31	0.05	0.3	0.2	23.6	0.56	0.63	0.2	0.344	0.404	21.4
			69.7	25.4	0.50	1.80	2.60	90.8	2.73	4.14	1.02	1.31
378/17B	6.14	500	5.26	1.72	0.07	0.3	0.2	25.5	0.115	0.104	0.02	0.2	0.404	20.8
			61.7	33.2	0.70	1.80	2.60	97.9	0.55	0.69	0.10	0.76
379/17B	6.18	510	3.24	1.64	0.05	0.1	0.2	19.4	0.08	0.104	0.01	0.112	0.296	15.2
			51.8	43.1	0.72	0.82	3.56	97.9	0.54	0.92	0.08	0.56

¹ Hard, hardness.

).

Low-TDS samples from Derbeke Depression (Region 2) are more acid, with pH varying from 6.07 to 6.18, than high-TDS samples from Regions 1 and 3. In the Dulgalakh River waters, pH increases to 6.72, and in the Dogdo River, pH rises to slightly alkaline at 7.71 due to the presence of carbonates (Ca²⁺ = 27 mg/L).

4.2. Stable Water Isotopes

High continentality explains a depleted isotopic signature in sampled water bodies, with mean values of δ¹⁸O = –20.6‰ ± 0.9‰, and δ²H = 149.7‰ ± 4.9‰. All but one sample plot was above the Global Meteoric Water Line (GMWL) (Figure 3 and

Table 4. Major ions content in samples from the Yana River headwaters (descriptions and location, see Table 2); numerator shows concentrations in mg/L, denominator, in %%-eq.

Lab #	δ18O, ‰	δ2H, ‰	dex, ‰	Lab #	δ18O, ‰	δ2H, ‰	dex, ‰
VKY17-001	−21.61	−154.9	18.0	VKY17-007	−20.01	−142.7	17.3
VKY17-002	−20.43	−147.5	15.9	VKY17-008	−20.91	−154.5	12.8
VKY17-003	−19.98	−146.4	13.4	VKY17-009	−18.94	−142.7	8.8
VKY17-004	−20.32	−149.6	13.0	VKY17-010	−19.86	−146.0	12.8
VKY17-005	−20.65	−151.1	14.2	VKY17-011	−21.94	−155.6	20.0
VKY17-006	−20.99	−148.5	19.4	VKY17-012	−21.87	−157.0	18.0

), suggesting important input from local precipitation. One sample, taken from a moist peat bog depression near the Lugovoye Lake shore (VKY17-009), exhibits evaporative enrichment effect with d-excess (d_ex = δ²H − 8 · δ¹⁸O) value d_ex = 8.8; in other samples, d-excess values vary between 12 and 20. The closest GNIP (Global Network of Isotopes in Precipitation) station is Yakutsk, Siberia, about 500 km southwest of the study area, so its information cannot be used for comparison in the scope of this study. Local meteoric water line cannot be established due to highly variable isotopic composition, inherited from mixing of meteoric precipitation with groundwaters of different origin and residence time (and icings).

4.3. Trace Elements & REEs

Variations in chemical composition follow regional pattern, but also the stream size (Figure 4a, Table S1). Samples cluster into three groups, uniting (1) the Dulgalakh and Dogdo Rivers, (2) lakes and streams from the Derbeke Depression, and (3) Narimchiki Lake, clustering separately from the latter. Rivers are clustered together due to the high content of Ca, K, Mg, S, Li, Sr, Mo and Sb, low Al, Y, Mn and Zr, but also Fe content below the detection limit. Lakes and streams from Region 2, inversely, show high concentrations of Be, Fe, V, Mn, Zn, Th, and low concentrations of Ca and S, notably. Narimchiki Lake sample 58B/17 plots as a separate cluster due to higher Na and lower Al and Si content, but a significant discrepancy in Na content between two analyses is noted: 50 µg/L by capillary ion electrophoresis, and 5000 µg/L by ICP-AES. Sample contamination cannot be excluded. However, it is unlikely since other elements in this sample exhibit reasonable behavior. Moreover, this sample differs from other samples taken in this region. A misprint in the datasheet provided by the analytical lab is also a plausible explanation.

The total REEs content is controlled by pH, with higher total REEs, from 0.8 to 1.1 µg/L in Region 2 samples with pH between 6.0 and 6.2, close to 0.5 µg/L in Region 1 sample with pH = 6.7, and is substantially lower at 0.06 µg/L in the Dogdo River sample, Region 3, with pH above 7.7.

The NASC-normalized REEs distribution in sampled waters is highly regional, showing three distinct patterns in relative abundance of light (LREEs) and heavy (HREEs) rare earth elements (Figure 5). The Dulgalakh River sample, Region 1, is depleted in HREEs, as shown by high La_N/Er_N = 2.1, La_N/Yb_N = 2.5, and shows unexpected Sm and Gd peaks: positive Sm anomaly Sm_N/Sm_N* = 7.3, and positive Gd anomaly Gd_N/Gd_N* = 4.0, calculated according to [86,87], respectively. Here and below, subscript N refers to NASC-normalized values, and * designates a ‘natural’ baseline REE concentration, against which the anomaly is calculated. All samples from Region 2 plot around the same pattern, with high HREEs content, La_N/Er_N and La_N/Yb_N below 0.4, low content and negative Ce anomaly: Ce_N/Ce_N* = 0.6…0.7, according to [87]. In common between these two patterns is a negative Tm anomaly [88], Tm_N/Tm_N* = 0.65…0.75 (Figure 5). In contrast, the Dogdo River sample, Region 3, shows significant depletion in LREEs, with La_N/Er_N = 0.32, La_N/Yb_N = 0.38. A pronounced negative Ce anomaly is found, Ce_N/Ce_N* = 0.4, and Eu content below detection limits, which could also imply its negative anomaly.

5. Discussion

5.1. Water Runoff Sources

The main objective of this study, besides the first in-depth description of water chemistry of the upper Yana River basin, is the multi-proxy assessment of potential water runoff sources in the three studied region. The dataset is limited in space and time, covering only late summer period, when the river-feeding sources are potentially most diverse. The possibility to perform a multi-proxy analysis based on three small datasets fully describing basic water chemistry was tempting, but the final result is of local importance.

5.1.1. Region 1: Dulgalakh and Nelgese River Basins

Rivers of the region plot above the GMWL and, with δ¹⁸O below −21.5‰, have the most depleted isotopic signature, and d-excess ranging from 18 to 20‰. On the δ¹⁸O vs. d-excess diagram (Figure 3b) these samples plot at the upper left corner, that can be indicative of sub-permafrost groundwater sources [19], or winter and spring precipitation [40]. In the sampling point, the Dulgalakh River had ultra-low TDS, about 50 mg/L, unsupportive of important sub-permafrost groundwater contribution, or evidencing its heavy dilution by precipitation or supra-permafrost groundwater. The absence of Ce anomaly (Ce/Ce* = 0.9) evidences short residence time in open channel, and points at runoff sources close to the sampling point. Sampling point location within the icing glade suggests that groundwater discharging within the glade, and potentially feeding the icing throughout winter, is from shallow supra-permafrost aquifer adjacent to the valley section. The presence of sulphates, high Nd/Yb = 19.9 and Ce/Mn = 159 ratios, point at relatively high filtration rates in an unconfined aquifer.

5.1.2. Region 2: Derbeke Depression

Region 2 samples have ultra-low TDS below 25 mg/L, a d-excess slightly above 10, and slightly acidic pH just above 6.0 suggest rainfall exhibiting minor modification at the rock interface as the dominant water source in sampled lakes and streams of the Derbeke Depression. Their chemical composition is resembling, since the samples represent two lakes, a tributary of one lake and the drainage from the other, but isotopic signature shows diverging patterns. Narimchiki Lake, its drainage stream and an inter-bog depression, Omchikandya River and a tributary to the Lugovoye Lake show similar isotopic features, d-excess between 12.8 and 14.2‰, presumably close to local precipitation, while the Lugovoye Lake itself and water from the soil pit close to the bank shore have higher d-excess, between 17.8 and 19.4‰, closer to fractionated groundwater. The REE distribution (Figure 5), convex-up with enrichment in MREEs, is indicative of the atmospheric water origin [89]. Water chemistry suggests migration pathways in the organic part of the soil profile; presumably, longer retention in the Narimchiki Lake catchment that is larger, and shorter pathways in the smaller Lugovoye lake catchment. Low content of highly soluble Ba and U, Ba/Th, and U/Th ratios below unity, and Nd/Yb ratio above 7 suggest fast subsurface water transfer and preferential REE binding to Fe(Al)-oxide colloids. Relatively high Fe content, 0.2 to 0.4 mg/L, is owing to leaching from organic soils and minor peat bogs. Hillslope water tracks are surficial runoff pathways extremely abundant in the region that enhance connectivity between rainfall and lakes and creeks.

5.1.3. Region 3: Dogdo River Basin

High pH, alkalinity and Ca²⁺ concentration suggest leaching of carbonates as a dominating process controlling water chemistry. Low total REE content and a marked depletion in LREEs are typical for carbonate rocks; besides, alkaline waters control REE concentration by complexation and co-precipitation with carbonate ions [90]. Increased Ba/Th and Ba/Ce ratios, as well as U/Th ratios exceeding 230 and Sr/Gd ratios over 29, point at deeper groundwater sources A pronounced Ce anomaly (Ce_N/Ce_N* = 0.4) implies longer residence times in both topsoil and open channel flow, and hence distant sources are postulated. However, it is also known to be a typical signature of carbonate dissolution [91]. A high U content, over 0.7 μg/L, suggests a potential input from calcic skarns of the middle reaches of the Dogdo River, but intrusive bodies all across the catchment could serve as alternative sources. Trace elements point to two major water source areas underlain by sandstones: the Dodgo River headwaters (high Mo and Rb), and the Sinekandya River, a left tributary of the lower Dogdo River (high As and Sb). Finally, suprapermafrost groundwater passing through fissured Middle and Late Jurassic sandstones is assumed to be a dominant water source in the Dogdo River.

5.2. Icing Discharge Signatures

Icings were observed in the Dulgalakh and Nel’gese River basins during field survey, and they are known to be present on the Dogdo River. However, whether it is possible to establish a certain ‘icing trail’ in sampled rivers, remains questionable. During freeze-up and icing development, the icing body is expected to retain isotopic signature of feeding groundwaters that can be highly variable. In central Yakutia and Lena-Viluy interfluve, a low d-excess in supra- and intra-permafrost springs is observed, and resulting icings have negative d-excess and significant evaporative enrichment [92]. As we discussed above, one sample in the active icing glade shows water depleted in ¹⁸O and with high d-excess, potentially highlighting the role of snow and low temperature freezing. Enrichment in Li, Mg and Sr is possible through cryogenic precipitation and remobilization upon thawing. The pattern of REEs interaction with carbonate precipitates is unclear; selective complexation and co-precipitation of HREEs with carbonates is expected, and upon carbonate dissolution upon meltdown, HREEs might be preferentially reintroduced to streams. However, relative importance of this process might be regional, related to carbonate rocks occurrence, and overridden by other processes

5.3. Sm and Gd Anomalies

Distribution of REEs shows several significant anomalies, of which some were not expected to occur in this remote location, namely high positive samarium and gadolinium anomalies. No published regional data exist on patterns in REEs content to be used as a baseline reference, but for both elements, their concentrations exceeded those observed in industrially-developed regions globally. The Sm content, 109 ng/L, is higher than in the heavily used and anthropogenically contaminated Rhine River, up to 31.2 ng/L [86]. The Gd content in natural waters is reported to vary from 1 to 4 ng/L [93], while it is consistently above 40 ng/L in most samples taken in this study. Several potential origins of these anomalies can be proposed: (1) atmospheric deposition; (2) anthropogenic contamination; and (3) local geochemical anomaly.

Atmospheric deposition can not be neither proved or ruled out, since our study lacks coverage to either support or reject this hypothesis. However, the Verkhoyansk Range is an effective orographic barrier for wet deposition, and hence significant contamination is unlikely to occur on its leeward slope.

Anthropogenic Sm contamination is shown to be related to chemical industry producing catalysts for petroleum refining [86], and Gd contamination is frequently associated with hospital effluents and MRI imaging contrasts [94]. Both chemical industry and well-equipped hospitals are absent in this remote and uninhabited area, and cannot cause the observed anomalies. Human presence in the area is limited to a small community, Seben-Kyol, with about 750 inhabitants, 70 km to the west, and Vertikalnoye silver-lead mine in an adjacent river basin, downstream the sampling point. Seben-Kyol operates a dirt runway for small Antonov-2 planes and helicopters, and a non-paved road connects it to the mine. Direct impact from the community is unlikely, since it is buffered from the sampling point by two large lakes, though this cannot be fully excluded—a network of winter roads is operational each year, and fuel leakage or other contamination episode will go unnoticed due to the remoteness and absent environmental control.

Chemical composition of local rocks is unknown, but according to the state geological map of mineral deposits [95], regional enrichment in Hg, Sb and Au is noted, but no major REEs accumulation. Other elements do not show particular anomalies in this sample; the presence of Cd, 0.05 µg/L, and elevated Li content, 1.4 µg/L, are shared with other regions. Cadmium is widely present in groundwater globally [96,97], and we assume the observed Cd, as well as associated REEs, to be geogenic, originating from rising sub-permafrost groundwater through an open talik. The sample was taken within an icing glade, supporting geogenic origin of Sm and Gd anomalies. In non-permafrost region, southwest England, groundwaters from Devonian metasedimentary rocks, greenshist grade, showed LREEs enrichment and Sm and Gd peaks [98]. Greenshists are described in the Verkhoyansk Range, in the vicinity of the sampling site [99], and could have affected water chemistry in the Dulgalakh River.

Although both elements are found in safe concentrations across the studied region, they are, however, significantly above the global background, and Gd has high bioaccumulation capacity. Future detailed samplings should confirm and delimit this anomaly. Fish tissue accumulation studies may be useful to understand and reduce potential population exposure to high Gd levels.

6. Conclusions

The studied region of the Verkhoyansk Range is remote and hard-to-access. Hence, regional water chemistry data are scarce. In total, 26 samples were collected during the helicopter-based reconnaissance survey. The first data on trace metals and rare earth elements in natural waters of the Upper Yana River are presented in the manuscript.

In late summer, the sampled water is of the HCO₃⁻—Ca²⁺—Mg⁺ type, and ultra-low TDS, from 20 to 100 mg/L. Multi-proxy analysis allowed for highlighting water sources in three distinct regions in the upper Yana River basin. At the northern slope of the Verkhoyansk range, the Dulgalakh and Nelgese Rivers, the dominant runoff pathway is a shallow subsurface flow in the unconfined aquifer. In the Elga Highland headwaters, the rainfall is delivered to lakes and streams via the shallow active layer and surficial flowpaths (hillslope water tracks). At the Chersky Ridge flank, the Dogdo River is fed by precipitation and groundwater percolating through fissured carbonates. Specific ‘icing trail’ in chemical composition can include high d-excess, above 18, and enrichment in Li, Ca and Mg, but icings should also be sampled as an important end-member in order to understand the evolution of groundwater chemistry upon discharge, freezing, and further meltdown. Interestingly, Sm and Gd anomalies were discovered in the Dulgalakh River; both elements are currently attracting attention as emerging contaminants. Greenshist dissolution represents potential geogenic source of both Sm and Gd, while fuel leaks, mining industry effects, and historical mining legacy are possible anthropogenic sources.