Black Sea hydroclimate and coupled hydrology was strongly controlled by high-latitude glacial climate dynamics

The Black Sea experienced pronounced millennial-scale changes in temperature and rainfall during the last glacial coinciding with Dansgaard-Oeschger cycles. However, little is known regarding the amount and sources of freshwater reaching this inland basin. Here, we present detailed ostracod δ18O data from the glacial Black Sea showing subdued Dansgaard-Oeschger cyclicity and four prominent longer-term saw-tooth shaped Bond-like cycles. We propose that the δ18Oostracods signature primarily reflects changes in the atmospheric circulation in response to the waxing and waning Eurasian Ice Sheet. The millennial-scale ice sheet variations likely resulted not only in latitudinal migrations of atmospheric frontal systems but also in shifts of dominant moisture sources for the Black Sea. Heavier isotopic precipitation arrived from the North Atlantic-Mediterranean realm during the warmer interstadials and lighter isotopic precipitation from the Eurasian continental interior during the colder stadials. The subdued Dansgaard-Oeschger variability likely reflects an integrated precipitation signal additionally affected by the long mixing times of the large Black Sea volume up to 1,500 years as suggested from hydrologic-isotope-balance modelling. Moisture sources to the Black Sea changed in response to atmospheric frontal displacements driven by Eurasian Ice Sheet dynamics during the last glacial period, according to analyses of ostracod oxygen and strontium isotope data from Black Sea sediments.

T o understand hemisphere-wide climate variability during glacials, profound knowledge of former oceanic and atmospheric coupling in a zonal and meridional perspective is essential requiring well-preserved palaeo-archives that are widely distributed over the oceans and continents. While especially records along the North Atlantic recording last glacial orbital-and millennial-scale climate variability are available, high-resolution continental records from the last glacial are still sparse. However, just the deeper continental palaeoclimate records may allow retracing high to mid-latitude teleconnection patterns of the atmospheric circulation system tightly coupled to the size of the Eurasian Ice Sheet during glacials. The former glacial isolated Black Sea, bridging the North Atlantic-Eurasian corridor, is one of such inland candidates. It could help to extend our view of high-to midlatitude atmospheric coupling during glacial climate variability.
During the relatively mild period within the last glacial between 57 and 29 ka BP, the Marine Isotope Stage 3 (MIS 3 1 ), temperatures fluctuated on a millennial time scale during the so-called Dansgaard-Oeschger cycles (DO cycles) with abrupt warming into interstadials and gradual cooling into stadials across the North Atlantic-Mediterranean corridor [2][3][4][5][6][7] . During MIS 3, the largely enclosed, lacustrine Black Sea experienced pronounced DO climate and environmental variability as well, which were tightly linked to changes in the North Atlantic Meridional Overturning Circulation (AMOC) and temperature variations in Greenland [7][8][9] . Temperature amplitudes between stadials and interstadials in the Black Sea were as high as 4°C 7 . During the warmer interstadials, the amount of rainfall and riverine freshwater supply into the basin increased 8,9 . The increased freshwater contribution resulted in increased lake level and lower salinity, whereas higher nutrient supply and warmer conditions favoured primary productivity. In Northern Anatolia, temperate forest expanded and wind strength was reduced 8,9 . During the colder stadials, the hydroclimate gradually changed from humid to arid conditions, which resulted in increased evaporation, a lowered lake level, and higher salinity. Along with reduced riverine supply of nutrients, primary productivity decreased. During stadial winters, temperatures were much lower than during interstadial winters resulting in ice formation at the southern coast of the Black Sea 'Lake' [7][8][9][10] . On land, the vegetation in Northern Anatolia became dominant in xerophytic steppe 8 . In a long-term perspective, similar changes in temperature and especially in rainfall amounts on an orbital timescale are superimposed to the DO-cycles with much more humid conditions during the first half of MIS 3 8,9 . During this time window it is supposed that the enhanced supply of freshwater by rainfall and rivers caused a rising lake level of the Black Sea that probably resulted into an overflow into the Marmara Sea until ca. 42 ka BP 9,11-13 . Previous studies have contributed to the knowledge about the changes in temperature and rainfall amount and the resulting environmental implications for the Black Sea 'Lake' and Northern Anatolia during MIS 3 [7][8][9][10] . In addition to the Black Sea sediment record, the δ 18 O records of speleothems from Northern Anatolia document these climate cycles exceptionally well and are suggested to reflect changes in the δ 18 O composition of the Black Sea surface water 14 . However, there is no clear explanation for the underlying processes and mechanisms of these apparent changes in the δ 18 O composition of the Black Sea surface waters during the intermediate climate state of Marine Isotope Stage 3 (57-29 ka BP). Therefore, it has to be unravelled how these changes relate to the variability of temperature, the amount and sources of rainfall and riverine supply as well as to the hydrological activity in the Black Sea itself. With a deeper knowledge of the coupling between hydroclimate variability and hydrological activity, it may help understanding continental high-to midlatitude teleconnection patterns of the atmospheric circulation during glacials. Here, we present a new continuous highresolution δ 18 O record obtained from benthic ostracod carbonate shells from a well-dated SE Black Sea sediment core in order to investigate the hydrological dynamics in the former Black Sea upper water column and their relation to migrating air mass trajectories subjected to hydroclimate reorganisations. In combination with a well-dated δ 18 O record from speleothems of the Northern Anatolian Sofular Cave that has been suggested to mainly reflect the δ 18 O composition of the former Black Sea surface water during the period of interest 14,15 , we will propose potential atmospheric-hydroclimatic processes explaining the similar and diverging trends in both records. We further will elaborate on why the δ 18 O record has a significant orbital timescale component and how this is related to the Eurasian Ice Sheet (EIS) dynamics and changes in atmospheric circulation patterns over Eurasia. A simple hydrologic-isotope-balance model will show that the δ 18 O ostracods pattern reflects an integrated central European precipitation signal with its imprint strongly affected by hydrological activity and mixing times of the large water volume of the Black Sea. A new 87 Sr/ 86 Sr ostracods record obtained from the same sediment core, will serve for estimating changes in the dominant source of riverine supply that might be partially associated with meltwater release from the EIS in turn potentially affecting the δ 18 O composition of the Black Sea.
Today, the Black Sea is a semi-enclosed, brackish basin and represents the largest anoxic water body on Earth 16,17 . It receives freshwater by rivers and rainfall and saltwater by the inflow of water from the Mediterranean and Marmara Seas through the Dardanelles and Bosporus straits 17 . Due to the higher input of freshwater by precipitation (P) and riverine runoff (R) compared to the loss of water by evaporation (E), the modern Black Sea has a positive water balance (P + R > E) 17 .
Controlled by the northern hemispheric atmospheric circulation over Europe, the climate in the Black Sea region is mainly controlled by the Polar Front Jet, the mid-to high latitude westerlies, the Subtropical Jet, and the Siberian and the mid-latitude subtropical high pressure systems 18,19 (Fig. 1). The dynamics and intensity of Mediterranean cyclones arriving in the Black Sea are affected by the East African and Indian monsoons in combination with the Siberian high-pressure system and relate to the positions of the Polar Front and Subtropical Jet 18,19 . With dry summers and wet winters, the Black Sea climate is generally under Mediterranean influence and its interannual variability is strongly linked to the North Atlantic Oscillation [20][21][22][23][24][25] . Mediterranean atmospheric moisture mainly originates from the North Atlantic, but evaporation of Mediterranean waters also strongly controls the isotopic composition of rainfall in the Black Sea drainage basin 26 .
Northeastern Anatolia and the eastern Black Sea are the most summer rain-laden Black Sea regions due to orographic precipitation in the Pontic and Caucasian Mountains [27][28][29] , whereas the Black Sea itself forms an additional moisture source for rainfall 14,15,24 .
On an over-regional scale, the δ 18 O composition of rainfall across Eurasia decreases progressively from the mid-latitudes to northern latitudes and generally eastwards to the continental interior ( Fig.1) due to the Rayleigh condensation [30][31][32] . The δ 18 O composition of rainfall in the Black Sea region ranges today from −4‰ to −8‰ 33 , whereas other studies reported even values of around −9‰ 34 . Additionally to the rainout effect with depletion towards northeast and the continent, the altitude of the NE Anatolian Mountains also affects the rainfall isotopic composition leading to values of about −10‰ 26 . Today, in the Black Sea the δ 18 O composition of the surface water ranges between −2.44‰ and −2.84‰ 15,33,35,36 , whereas the deeper water is slightly heavier and reveals a composition between −1.65‰ and −1.8‰ 33,35,36 . However, it has to be noted that during strong precipitation periods, the oxygen isotopic composition of the surface water can decrease to −7.1‰ 33 , demonstrating the pronounced sensitivity of the Black Sea to variations in freshwater input. With 50-70% of the total riverine input, the Danube river is the main freshwater source from all rivers draining into the Black Sea today 17,37 . Near the river mouth, the Danube reveals a δ 18 O composition of −10.5‰ 36 . While the freshwater contributed by rivers and rainfall is isotopically highly depleted, the high-salinity Mediterranean water entering the Black Sea has a δ 18 O composition of about +1.8‰ 33 .

Results
The presented data derive from the well-dated 10      . Apart from some outliers during the middle part of MIS 3 between 40 and 45 ka BP, the 87 Sr/ 86 Sr ostracods pattern bears resemblance to the long-term pattern of the Sr/Ca ostracods record, with its minimum level interpreted to reflect a freshening phase during MIS 3 9 (Fig. 3g).

Discussion
Potential factors changing the Black Sea δ 18 O ostracods composition during the last glacial. The high-resolution δ 18 O ostracods record presented here on a centennial timescale reveals at first glance a subdued nature of DO variability suggesting that it neither reflects simply temperature nor rainfall amount since both abruptly changed during last glacial stadial-interstadial cycles in the Black Sea region [7][8][9] . The variations in the δ 18 O ostracods record represent a combined response to changes in the δ 18 O composition of the water and in the ambient temperature during the growth of the ostracod shells. The δ 18 O composition of the glacial lacustrine Black Sea water column, in turn, was affected by changes in a) evaporation of the lake surface and b) the amount, sources, and δ 18 O composition of freshwater input from direct precipitation at the lake surface and indirectly via rivers arriving from the Black Sea drainage basins. The latter isotopic signature comes from a mixture of precipitation and potential meltwater from adjacent glaciers and ice sheets feeding the rivers 38  When affected only by temperature, the range of about 3.4‰ in the δ 18 O ostracods would suggest temperature amplitudes of about 13.6°C during the investigated period 43,44 . However, temperature changes between stadials and interstadials reached only up to 4°C during the period 64-20 ka BP as revealed by reconstructions of TEX 86 -based annual lake surface temperatures 7 (Fig. 3d). Similarly, estimations of temperature changes in the intermediate water depth of~400 m by Mg/Ca ostracods imply temperature amplitudes of 2-4°C during the DO-cycles 7 . This circumstance and because the record a) lacks clear millennial-scale DO variability and b) shows generally higher values during warmer periods although the opposite pattern (negative correlation) would be expected, the temperature impact on the δ 18 O ostracods was most likely very low 43,44 .
By considering a relation of 0.25‰ decrease in δ 18 O per rise of 1°C in temperature 43,44 , we corrected the δ 18 O ostracods values by removing this temperature effect using the TEX 86 -based lake surface temperatures from the same sediment core 7 to unfold a potential DO variability not coupled to temperature. Although this lake surface temperature correction relative to the conditions at the ostracod's habitat depth at 418 m represents an overcorrection, the temperature-corrected δ 18 O ostracods values are very much the same as the non-corrected values (Fig. 3c) confirming another superordinate environmental factor for the δ 18 O ostracods composition of the Black Sea 'Lake'.
Similar to an enhanced temperature, increased input of freshwater by rivers including meltwater discharge would result in lowered δ 18 O values. Previous studies have demonstrated a pronounced variability in the amount of riverine input on orbital-and millennial timescales from different proxies from the same sediment core presented in this study 8,9 . Since riverine discharge includes both direct rainfall and meltwater, we try to unravel the impact of these different freshwater sources on the δ 18 O ostracods record in this chapter. Irrespective of the freshwater source, an increase must have lowered the δ 18 O composition of the water in the Black Sea basin. On a millennial timescale, rainfall was much higher during the warmer interstadials as indicated by an increased arboreal pollen abundance reflecting the expansion of moisture-dependent temperate forest vegetation 8 . The related riverine discharge likely lowered the surface water salinity as mirrored by decreasing Sr/Ca ostracods and increased freshwater dinoflagellate assemblages during interstadials 8,9 (Fig. 3g). Furthermore, the amount of leaf wax derived odd-numbered long-chain n-alkanes along with the sedimentary ln(K/Ti) ratio, both reflecting relative changes between riverine and aeolian input, support the observation of an amplified riverine discharge during interstadials 9 .
In addition to the DO-patterns, the long-term patterns of the aforementioned proxies suggest a pronounced orbital-scale component with increased freshwater supply during the first and warmer half of MIS 3 strongly resembling the pattern of the coevally reduced Eurasian Ice Sheet volume 45 (Fig. 3i). Therefore, in addition to rainfall, increased meltwater discharge originating from the Eurasian Ice Sheet could have affected the hydrology and salinity of the former Black Sea 'Lake'. Similar as during other meltwater pulses during the last two glacial terminations and the deeper penultimate glacial, the Dnieper could have played an increasing role in freshwater discharge relative to the Danube 39,46,47 . A reliable tool for estimating changes in river sources for the glacial Black Sea is the Srisotopic composition of ostracods 39,47,48 , because it depends on the catchment geology in the respective drainage regions 48,49 . Our new 87 Sr/ 86 Sr ostracods record (Fig. 3h) shows generally less radiogenic values when compared to the last two glacial-interglacial terminations (T I 48 : 0.709102; T II 47 : 0.709452) and the penultimate glacial (MIS 6 39 ; max. 0.709405) and lowest values are reached during MIS 3 (0.7085460) with a slightly elevated glacial background during the colder periods (MIS 4 and 2: 0.708789). These values indicate an increased contribution of Dnieper and Don relative to Danube 50 including meltwater from the partly disintegrating Eurasian Ice during the earlier part of MIS 3 45 . When disregarding runoff from the remaining smaller rivers, our estimations imply that the contribution by Dnieper relative to Danube could have increased up to 71% during MIS 3, whereas during MIS 4 and 2 the relative contribution by Danube would have amounted up to 78% (Fig. 3h). A certain contribution from Caspian waters ( 87 Sr/ 86 Sr: 0.70821-0.7085) via spill out along the Manych depression (Fig. 1) is possible as well [50][51][52] . Although not explicitly suggested from the Sofular Cave speleothem δ 18 O record, the dinoflagellate assemblages from the Black Sea sediment core would be in line with a potential Caspian connection 8,15 .
Despite the clear evidence for early MIS 3 meltwater contribution, the amount released by the medium sized MIS 3 ice sheet 45 was obviously insufficient for decreasing the δ 18 O composition of the Black Sea water column as it was the case for the meltwater pulses of the last two glacial terminations and the deeper penultimate glacial 38,39,47 . Comparable trends in the δ 18 O ostracods , 87 Sr/ 86 Sr ostracods , and EIS volume records (Fig. 3b, h, i), however, suggest changes in the atmospheric circulation over the North Atlantic and Eurasia as a potential trigger of the observed patterns. Alike the Black Sea ostracod record, the δ 18 O values from the speleothems of Sofular Cave in north-western Anatolia are also heavier during the more humid interstadials than during the stadial periods 14,15 (Fig. 3e). While the speleothem data are suggested to be closely linked to the surface water conditions of the former Black Sea 'Lake' 14 , it is conceivable that the δ 18 O Sofular Cave was affected not only by the Black Sea as moisture source 14 , but incorporated a variable portion of heavier atmospheric moisture signal originating in the (eastern) Mediterranean Sea and the North Atlantic.
Since a retreated Eurasian Ice Sheet was most likely associated with a northward migration of the atmospheric polar front, precipitation patterns and moisture sources may have changed as well. The northward migration of the frontal systems over Eurasia possibly caused a progressively stronger influence of the North Atlantic/Mediterranean moisture sources to the precipitation in the Black Sea region and could explain a shift towards a heavier δ 18 O signature. Arppe and Karhu 53 pointed out that present precipitation patterns over Europe with a westerly moisture source and a progressive isotopic depletion in the heavier isotopes of air masses towards northeast existed also during MIS 3. Based on δ 18 O analyses of mammoth tooth enamel and palaeogroundwaters across Europe, the authors reconstructed the oxygen isotopic composition of rainfall for the period 52-24 ka BP, showing ca. 2‰ lighter values relative to present 53 . Although there are no reconstructions of δ 18 O-precipitation patterns over Eurasia available for more discrete intervals or even on a millennial timescale, the general reconstruction for the MIS 3 implies lighter δ 18 O values in SE Europe when the Eurasian Ice Sheet was present (MIS 3) than it was absent (today) 53 . This supports the idea of the MIS 3 Eurasian Ice Sheet dynamics as an important forcing factor in modulating the Eurasian atmospheric circulation pattern and moisture sources of the Black Sea region.
According to reconstructions by Bintanja and van de Wal 45 , the Eurasian Ice Sheet began to expand again since ca. 40 ka BP (Fig. 3i). Accordingly, a repeated southward migration of the atmospheric polar front is indicated from 40 ka BP by increasing wind intensity and drier conditions in SE Europe as reconstructed from grain-size distribution of loess-palaeosol sequences in the Lower Danube Basin 54 . The authors argue for an increasing role of the Siberian/Eurasian high-pressure system in modulating the continentality towards the peak glacial phase (Fig. 1). Such a gradual change likely reduced the North Atlantic/Mediterranean sourced precipitation in the Black Sea region causing a strong negative trend in the δ 18 O ostracods record due to increased contribution of rainfall sourced in the Eurasian continental interior. On the millennial scale, the stronger imprint of Heinrich stadials compared to the regular stadials in the δ 18 O ostracods record would then similarly point to a major reorganisation of the atmospheric circulation with a stronger continentality and much lighter δ 18 O ostracods . This goes in hand with the expected direct effect from the condensation temperature of moisture sources above the cave and more generally the Black Sea drainage basin with lighter (heavier) δ 18 O at low (high) stadial (interstadial) temperatures.
Changes in the Black Sea 'Lake' hydrology during the last glacial. On the millennial scale the δ 18 O ostracods record is not following simply the DO pattern, but appears to be subdued and slightly delayed regarding the peak interstadials. This may be caused by variable mixing times of the relatively large Black Sea 'Lake' water body in response to the hydrological balance of the basin. Similar to our MIS 3 record, δ 18 O ostracods records from the western Black Sea of the last glacial maximum and Termination I show rather smoothly increasing values and lack a clear shift during the Bölling-Alleröd warming and the cold Younger Dryas stadial 38,55,56 .
The present Black Sea water residence time is in an order of 1000 years 57 but it might have varied significantly during the last glacial period. It can be assumed that the mixing time in a basin of that size and geometry could have lasted as long as or even longer than an average DO cycle resulting in a significant modification of an original freshwater input signal. The mixing time during the late MIS 3/MIS 2 increased potentially to several millennia in response to the generally aridification trend towards the LGM with a negative hydrological balance in the basin. Higher evaporation and a decreasing Black Sea 'Lake' level is indicated e.g., by increased salinities as reconstructed by Sr/Ca ostracods 9 , sequence-stratigraphic evidence of an LGM lake level low stand in the Danube fan area 58 , and significantly higher reservoir ages of up to 1450 years towards the LGM 10,59,60 .
Following the approach by Bahr et al. 38 performing a simplified hydrologic-isotope-balance model (HIBAL 61 ) for the Black Sea since the LGM, we simulated the δ 18 O composition of the Black Sea for the period 64-20 ka BP by considering potential changes in the amount of precipitation, runoff, and evaporation, which consequently involves changes in the mixing time of the Black Sea water. As a high-latitude rainfall source we used the δ 18 O NGRIP 62-65 record (Supplementary Methods, Supplementary Table S1, Supplementary Fig. S1). Supplementary Table S1 shows under which specific but simplified boundary conditions during four time slices our simulated δ 18 O record best mimic the measured δ 18 O ostracods record (Fig. 3b, c; Supplementary Fig. S1).
Considering a reservoir age of 1450 years for the LGM around 23 ka BP 10 lying in our fourth time slice 20-28 ka BP of our HIBAL estimations (Supplementary Table S1), mixing times for the other three times slices can be roughly estimated from the hydrologic input volumes and are~1500 years for 65-60 ka BP, 700 years for 60-45 ka BP, and~1000 years for 45-28 ka BP. These simulations support not only a high-latitude rainfall source, but also long and varying mixing times of the Black Sea water explaining the subdued nature of DO-variability in the intermediate water depth of the glacial, lacustrine Black Sea. For instance, during the earlier part of MIS 3 predominantly more humid conditions and a certain contribution by meltwater 8,9 suggest a generally more enhanced hydrological activity. These conditions caused a rising lake level that probably culminated into an overflow into the Marmara Sea 9,11-13 associated with reduced mixing times and lower reservoir ages 10,60 being half as high as during the ending MIS 4 and towards the LGM.
Interestingly, the Bond-cycle like pattern in the δ 18 O ostracods with a subdued DO-variability is remarkable similar to the pattern of global sea-level changes 66,67 (Fig. 3 b, c, f). It was proposed that gradual sea-level rises were related to interstadial warming and resulted mainly from melting of northern hemisphere ice sheets 66 . Obviously, the oxygen isotopic composition in the isolated Black Sea 'Lake' was not directly affected by such changes, but centennial scale changes in the dynamics of the Eurasian Ice Sheet and modified atmospheric teleconnections between the North Atlantic and Eurasia 7 might have had an additional imprint on the δ 18 O ostracods .
The new high-resolution δ 18 O ostracods record covering the last glacial from 64 to 20 ka BP along with the hydrologicisotope-balance estimations demonstrates that the hydroclimate regime in the Black Sea region was strongly coupled to the northern hemisphere climate variations. During periods of a retreated Eurasian Ice Sheet, the δ 18 O ostracods generally increase, which is likely related to the northward migrations of atmospheric frontal systems and a reduced continentality with an increasing influence from North Atlantic/Mediterranean moisture sources for the Black Sea region. Although the decreased 87 Sr/ 86 Sr ostracods values during the milder MIS 3 point to an increased meltwater contribution from the Dnieper, Don, and possibly Volga rivers, the effect is not evident in the δ 18 O ostracods record. The subdued nature of the Dansgaard-Oeschger pattern recorded in the δ 18 O ostracods we attribute to the long mixing times of the former Black Sea 'Lake' reacting sensitive to hydrological changes in the basin. We suggest generally shorter mixing times in an order of~700 years in an overflowing and fresher Black Sea in the first warmer and more humid part of MIS 3 and increasing mixing times up to~1500 years during the later part of MIS 3 and during MIS 2 and late MIS 4 with generally more arid conditions and a lowered lake level.

Methods
Sediment core and chronology. This study bases on the gravity core 25GC-1 recovered during RV Meteor cruise M72/5 in the SE Black Sea in 2007 ( Fig. 1; Archangelsky Ridge, 42°06.20′N, 36°37.40′E, 418 m water depth). The sediment core consists mainly of detrital mud intercalated by carbonate-rich intervals. The core chronology 10 (Fig. 2) bases on AMS radiocarbon dating, the Campanian ignimbrite tephra (Y5; 39.28 ± 0.11 ka 68 ), the Laschamp geomagnetic excursion (40.70 ± 0.95 ka BP 69 ), and proxy-tuning to the δ 18 O record from the North Greenland Ice core Project on the GICC05 timescale 62-65 (NGRIP). A study on cryptotephra horizons from the same core largely confirmed the NGRIP-tuned stratigraphy 70 . The sequence of interest in the present study covering the last glacial from the late MIS 4 to the MIS 2 (64-20 ka BP) derives from core depth 307-952 cm.
Ostracod isotope geochemistry. For the analyses of the δ 18 O and 87 Sr/ 86 Sr isotope signatures of the ostracods, well-preserved valves of adult specimens of Candona spp. were hand-picked from the wet-sieved and dried sediment fraction >150 µm.
For the oxygen isotopic composition of the ostracods, the carbonate valves from 447 samples in total were cleaned with a thin brush and deionised water under a binocular microscope. Then, the pre-cleaned ostracods (20-80 µg) were measured with a Finnigan MAT 253 isotope ratio mass spectrometer coupled with an automated KIEL IV carbonate preparation device (103% H 3 PO 4 at 70°C) at the GFZ Potsdam, Germany. Oxygen isotope values are given in delta notation relative to Vienna Peedee Belemnite. Repeated measurements of the international reference material NBS 19 and a laboratory internal standard (C1) assured an analytical precision better than ±0.07‰ for δ 18 O ostracods .
For the analyses of the strontium isotopic composition, 35 samples evenly distributed along the core with each 10-15 ostracod valves were cleaned with a few drops of H 2 O 2 (5%), methanol, and deionised water 48,71,72 . After dissolution of the valves in 300 µL 2 vol% HNO 3 , the samples were dried on a hot plate at 80°C. Strontium was separated from the sample matrix by using the HNO 3 -H 2 O technique in microcolumns filled with Eichrom © Sr-spec resin. Isotope analyses were performed with samples re-dissolved in 0.3 N HNO 3 and all isotopic beams measured simultaneously in static mode on a MC-ICP-MS (Neptune Plus, Thermo) at the University of Tübingen, Germany. Instrumental mass bias was corrected using exponential law and an 88 Sr/ 86 Sr ratio of 8.375209. External reproducibility for NBS SRM 987 was 0.710253 ± 0.000015 (2 SD, n = 16) for the 87 Sr/ 86 Sr ratio and total procedural blank was <95.8 pg Sr.

Data availability
The data used in the present study are compiled in the Supplementary Data 1 and are freely available online at the data repository Zenodo (https://doi.org/10.5281/ zenodo.4545579).