Exceptional warmth and climate instability occurred in the European Alps during the Last Interglacial period

Warmer temperatures than today, over a period spanning millennia, most recently occurred in the Last Interglacial period, about 129,000 to 116,000 years ago. Yet, the timing and magnitude of warmth during this time interval are uncertain. Here we present a reconstruction of temperatures in the Swiss Alps over the full duration of the Last Interglacial period based on hydrogen isotopes from fluid inclusions in precisely dated speleothems. We find that temperatures were up to 4.3 °C warmer during the Last Interglacial period than in our present-day reference period 1971 to 1990. Climate instability, including an abrupt cooling event about 125,500 years ago, interrupted this thermal optimum but temperatures remained up to 2.0 °C warmer than the present day. We suggest that higher-elevation areas may be more susceptible to warming relative to lowland areas, and that this may hold also for a future climate forced by increasing levels of greenhouse gases. Air temperatures in the European Alps during the Last Interglacial were up to 4 °C higher than at present, which could suggest amplified warming at high elevations, according to fluid inclusion data from Swiss speleothems.

T he Last Interglacial period (LIG), roughly equivalent to Marine Isotope Stage (MIS) 5e,~129-116 thousand years ago (ka), is currently the focus of increased interest as it provides a potential test bed for warmer conditions on Earth. In central Europe, previous work has found that summer temperatures may have been about 1-2°C warmer than present during the optimum [1][2][3][4] , with alpine regions possibly experiencing warmer temperatures of +4.3 ± 1.6°C 5 . Additionally, it is becoming more widely accepted that millennial-to centennial-scale swings in climate were intrinsic to this time period 6,7 . However, the magnitude, as well as the timing and spatial variability of climate during this period are still widely debated. This is largely due to the paucity of independently dated records that cover the entire LIG, with even fewer records having the necessary resolution to capture subtle climate oscillations. Because of this, quantitative temperature reconstructions are rare. The few available temperature reconstructions are from biogenic archives 6,8,9 , leading to potential biases. Further, a critical lack of age control for biogenic proxies has led to reported temperature fluctuations being poorly constrained in time. To address this, we reconstruct temperatures using fluid inclusions extracted from precisely dated speleothems. Fluid inclusions (fi), which are natural repositories of precipitation falling above the cave, are reliable indicators of paleotemperature [10][11][12][13] . The measurement of hydrogen isotopes (δD fi ) in speleothem fluid inclusion water can be performed with high precision 14,15 , allowing robust paleotemperatures to be calculated. Here, we provide a centennial to millennial-scale temperature reconstruction (dubbed SKR-FIT for Schrattenkarst fluid inclusion temperature) using δD fi fluid inclusion data from speleothems collected from two subalpine caves in Switzerland (Fig. 1). The SKR-FIT record indicates temperatures up 4.3 ± 1.4°C higher than modern-day (1971)(1972)(1973)(1974)(1975)(1976)(1977)(1978)(1979)(1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990) during the LIG and a 1.9 ± 1.4°C decrease over this thermal maximum.
Site setting. The caves are located on Schrattenkarst in the Melchsee-Frutt region of central Switzerland (46°47′N, 8°16′E) ( Fig. 1 and Supplementary Fig. 1), a karst region characterized by a high density of caves with a cumulative length of 65 km. The caves developed in different paleo-phreatic levels between 1300 and 2450 m above sea level (a.s.l.) 16 . The host rocks are limestones (Quinten Formation) and limestone-marl successions (Mürtschen Member of the Schilt Formation 17 ) of uppermost Jurassic to lowermost Cretaceous age. Speleothems from two caves close in proximity, Neotektonik Cave and Schratten Cave, were studied. The former opens at 1727 m a.s.l. and comprises 2.1 km of galleries; the latter is much larger (19.7 km) and its entrance is located at 1716 m a.s.l 18 . Three stalagmites (M37-1-16A, M37-1-16C, M37-1-23A), found broken and inactive, were collected 75-103 m behind and 27 m below the entrance of Neotektonik Cave where the rock overburden is 30-35 m. One stalagmite (M6-73-3) was found in Schratten Cave 475 m behind and 136 m below the entrance where the rock overburden reaches 54 m. The catchment of the cave is located at the transition from the subalpine zone with a few scattered stands of spruce to alpine grasslands. The natural tree-line in this area is located at about 1900 m a.s.l., but mostly lower due to century-long pasturing. The catchment is partly covered by thin, brown soil and partly barren karstified rock. The mean annual air temperature at the elevation of the cave is 2.9°C based on data from nearby meteorological stations (1990-2018; see "Methods" section). Mean annual precipitation at the study site is~2200 mm.

Discussion
The outstanding features of the SKR-FIT reconstruction are temperatures up to 4.3 ± 1.4°C higher than modern-day (1971)(1972)(1973)(1974)(1975)(1976)(1977)(1978)(1979)(1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990) during the LIG and a 1.9 ± 1.4°C decrease over this thermal maximum (Figs. 2 and 3). The warm temperatures corresponding to the LIG thermal maximum occur during two separate time periods, between 127.3 ± 0.7 and 125.9 ± 0.5 ka, and between 124.6 ± 1.0 and 124.1 ± 0.5 ka (Figs. 2 and 3). During these periods, the SKR-FIT indicates median temperatures 3.6 ± 1.4°C higher than modernday, with temperatures up to 4.3 ± 1.4°C higher (Figs. 2 and 3). This is significantly higher than lower elevation (<1000 m a.s.l.) pollenbased temperature reconstructions of 0-2°C above modern-day 2,19 . This contrasting temperature difference between pollen-based reconstructions and the SKR-FIT could be due to pollen records corresponding to either summer (July) or winter (January) temperatures 2 , although pollen records that record mean annual temperatures 19 also show similar low-temperature anomalies during the early LIG. Differentiating climate estimations between records could also account for the temperature difference, with pollen-inferred temperature reconstructions relying on the use of vegetation-climate transfer functions of key indicator taxon or statistical models linking climate variables and vegetation distribution 20 . Additionally, it is always of concern if pollen-inferred temperature fluctuations during the LIG are reliably constrained in time, which rely solely on biostratigraphy without radiometric age control. If the pollen records provide reliable paleotemperature estimates, then the difference in temperature may instead be a result of elevation, rather than proxy biases. While no high-elevation pollen-based temperature reconstructions exist in Europe for the LIG to directly compare a possible elevation dependence between records, the high-elevation Cesare Battisti Cave in the northern Italian Alps 5 (1930 m a.s.l., Fig. 1) agrees within uncertainty with the SKR-FIT record. Temperature reconstructions from this cave record an increase of 4.3 ± 1.6°C relative to 1961-1990, with the high-temperature increase, interpreted as a result of elevation-dependent warming, meaning high-altitude sites react more sensitively to warming than lowelevation sites 5 . The possible significant increase in temperature during the LIG in the European Alps compared to modern-day is consistent with Arctic regions, which show a similarly pronounced warming of 4-5°C during the LIG 21 (Fig. 3). Additionally, deposition of ice-rafted debris (IRD) along the margins of both east and southwest Greenland indicates that iceberg discharge persisted during the early LIG, and was possibly associated with meltwater events 22,24 ( Figs. 1 and 3). Further, a prominent red-colored IRD layer, remarkably similar to a red-colored IRD layer deposited during the 8.2 ka BP event, was identified in the Hudson Strait and is associated with the final collapse, and outburst flood, of the Laurentide Ice Sheet 25 . The red layer in MD03-2664 was assigned an age of 124.2-124.7 ka 23 . Terrestrial evidence of a cold event associated with C27 is ambiguous. The most convincing evidence of cooling was reported from Sokli in northern Finland, where pollen extracted from a sedimentary sequence indicate summer cooling of 2°C, which was linked to meltwater pulses into the North Atlantic 6 (Fig. 1). However, this cooling event is largely unconstrained in time, with the event possibly being related to cold event C28, a hypothesized meltwater event during the early LIG dated tõ 128.5 ka 7 . A cooling signal during the early LIG was also identified in a speleothem record from Baradla Cave in Hungary between 127 and 123 ka 26 (Figs. 1 and 3). However, due to low sample resolution in this interval, correlation to a specific event is unconstrained. Therefore, the low temperatures between 125.8 ± 0.5 and 124.6 ± 1.0 ka in the SKR-FIT record are the first well-constrained evidence of a cooling associated with cold event C27 in central Europe. This implies that the impact of C27 was much more widespread than previously thought, with the event possibly being comparable in magnitude to the 8.2 ka BP event 23 , which resulted in a 1.7°C mean annual air temperature decrease in central Europe 27 . After 124.1 ± 0.5 ka, the SKR-FIT record suggests gradually decreasing temperatures for the remainder of the record (Figs. 2 and 3). However, temperatures remained 0.5 to 1.8 ± 0.8°C above modern-day temperatures, even during the glacial inception between 122 and 115 ka 28 . Warm temperatures towards the end of the LIG are consistent with pollen and SST data from core MD95-2042 off southern Portugal which show persistent Eurosiberian and Mediterranean plant species, as well as warm SST's well into MIS 5d 29,30 (Fig. 1). It has been proposed, based on proxy and modeling studies [31][32][33][34][35] , that strengthened NADW ventilation during late MIS 5e resulted in continued tropical heat transport into the North Atlantic. While the higher-than-modern temperatures revealed by our record support the interpretation of continued heat transport into the North Atlantic during late MIS 5e, cooling conditions between 124.2 ± 0.5 and 115.8 ± 0.4 ka reveal a gradually deteriorating climate.
Our new quantitative data document exceptionally warm conditions during the LIG in the Alps, as well as a prominent climate reversal, when atmospheric greenhouse gas concentrations were comparable to pre-industrial values of the late Holocene. Although these two warm periods cannot be directly compared because of different orbital settings, the high temperatures and the mounting evidence of unstable nature of the LIG are alarming in the context of future anthropogenic greenhouse projections. In particular, our results suggest that higher elevation sites may amplify warming compared to lower elevation sites. This may imply accelerated warming at higher elevations with increasing greenhouse gases. Additionally, The SKR-FIT record provides compelling evidence that increasing warmth will likely be associated with future climate instability in Europe, which may be linked to a reduction in NADW formation from accelerated melting of the Greenland Ice Sheet.
U-Th ages. A total of 38 powdered calcite samples were manually drilled for U-Th dating under a laminar flow hood; 24 from M37-1-16A, 5 from M37-1-16C, 4 from M37-1-23A, and 5 from M6-73-3 ( Supplementary Fig. 2 and Supplementary Data 1). U-Th samples were processed at the University of Minnesota Trace Metal Isotope Geochemistry Lab and analyzed using a ThermoFisher Neptune Plus multi-collector inductively coupled plasma mass spectrometer equipped with an Aridus desolvation nebulizer, following the method of Shen et al. 36,37 . Ages are reported with 2σ errors in years before present, i.e., before the year 1950 A.D. A time-depth model was created in OxCal 4.3 using the Bayesian approach 38,39 . The duration of the hiatuses at 78 and 225 mm in M37-1-16A, and at 45 mm in M37-1-16C, was estimated using the time-depth model ( Supplementary Fig. 3).
Stable isotopes. A total of 2905 stable isotope locations were drilled using a Merchantek micromill (Supplementary Data 3). In M37-1-16A, samples were drilled every 0.25 mm from 0 to 78 mm, and every 0.5 mm from 78 to 884 mm, with a temporal resolution of~5-10 years. Samples in M37-1-16C were drilled every 0.3 mm, with a resolution of~5 years. In M37-1-23A, samples were drilled every 1 mm, with a resolution of~70 years. Lastly, in sample M6-73-3, samples were drilled every 0.5 mm, with a resolution of~14 years. Stable isotope samples were analyzed at the University of Innsbruck using a ThermoFisher Delta V isotope ratio mass spectrometer equipped with a Gasbench II 40 . Stable isotopes are reported in per mil relative to Vienna Peedee Belemnite (VPDB). Long-term analytical precision is better than 0.08‰ for both δ 13 C and δ 18 O (1σ).
Fluid inclusions. Speleothem fluid inclusion water isotopes were analyzed at the University of Innsbruck using continuous-flow analysis of water via hightemperature reduction on glassy carbon 15 . δD fi isotope ratios are given in per mil (‰) using the standard delta notation and are reported relative to the Vienna Standard Mean Ocean Water (VSMOW). We extracted 28 calcite blocks from 16  Fig. 2 and Supplementary Data 2), weighing between 2.5 and 3.5 g, from the central growth axis of stalagmites M37-1-16A, M37-1-16C, M6-73-3, and M37-1-23A. Replicates were produced at 8 of the 16 depths (Supplementary Data 2). See Dublyansky and Spötl (2009) 15 for details on the crushing procedure. The precision of replicate measurements of our in-house calcite standard is typically 1.5‰ for δD fi for water amounts between 0.1 and 1 μL. Because crushing of our calcite samples released up to 0.86 μL of water (mean 0.25 μL), the precision of 1.5‰ for δD fi was found to be adequate for this study. Temporal resolution varied between several hundred to several thousand years. The paleotemperature record of stalagmites M37-1-16A, M37-1-16C, M6-73-3, and M37-1-23A was reconstructed based on the modern-day regional water isotope-temperature relationship 41 . We assumed that the lapse rate was the same for the LIG, and used δ 18 O values of modern precipitation and corresponding temperatures from the Global Network for Isotopes in Precipitation (GNIP) station at Grimsel Pass. This station was selected as it is the closest station (~24 km), with similar elevation (253 m elevation difference), to Schrattenkarst. We opted not to incorporate paleotemperatures based on equations from 5,42 , as they show a low correlation between calcite δ 18 O and the SKR-FIT record (Supplementary Fig. 4). We would expect some degree of correlation between δ 18 O and the SKR-FIT record if we assume fractionation effects are minor, given that they are both derived from precipitation. Furthermore, the reconstruction based on the modern-day regional water isotope-temperature relationship shows temperatures above 0°C for the duration of the record, necessary for speleothem deposition. Equations from refs. 5,42 , on the other hand, show temperatures well below 0°C, which are unrealistic for speleothem deposition. Only δD fi values were used for calculating paleotemperatures for the following reasons: post-depositional processes can alter the original δ 18 O fi in fluid inclusion water and thus limit the use of δ 18 O fi for paleotemperature calculations 43 . And, δD fi is not affected by isotopic fractionation during calcite precipitation and remains unaltered as there is no hydrogen source once the water is entrapped in the calcite matrix. To determine δD for Grimsel Pass, we obtained δ 18 O measurements from the GNIP dataset and converted to δD, as the GNIP dataset is missing several years of δD measurements. The δ 18 O/annual air temperature slope was determined from the Grimsel Pass GNIP station for the period spanning 1971-1990. This time period was selected in order to limit Anthropocene warming effects, which generally rises sharply after 1990. Because it is unclear which transfer function might be appropriate for the LIG, we evaluate a range of possible transfer function scenarios for the LIG, between 0.6 and 0.7‰/°C, which represents the error of the Grimsel Pass transfer function. Next, a δD/T transfer function (5.2‰/°C for Grimsel Pass) was determined by converting the annual air temperature slope by a factor of eight. We used the global meteoric water line relationship (δD = 8δ 18 O + 10) to convert δD fi to δ 18 Supplementary Fig. 5). The mean annual air temperature (MAAT) at Schrattenkarst between 1971 and 1990 is 2.1°C based on linear regression of the regional lapse rate of 0.51°C/100 m from Grimsel Pass (Supplementary Fig. 6). Further, δD was corrected for the ice volume effect of −0.064‰ per meter of sea level rise 37 , with sea level during the LIG approximated using the mean global sea level reconstructed from core KL11 from the Red Sea 44 . The modern uplift rates in the study area arẽ 0.5 mm/yr 45 , which translates into~60 m of uplift since the LIG. Rock uplift rates due to the melting of the Last Glacial Maximum ice cap and concomitant erosion contribute significantly to the modern vertical displacement rates in the Alps 46 . We, therefore, regard the~60 m of uplift since the LGM as a likely maximum value. Based on the regional lapse rate, this equates to a temperature decrease of~0.3°C since the LIG. To account for this uplift component, we decreased the SKR-FIT values by 0.3°C (Supplementary Data 2). We do not account for changes in the seasonality of precipitation nor changes in atmospheric circulation patterns across Termination II and into the LIG. Fluid inclusion data were also obtained from a Holocene stalagmite from Neotektonik Cave (M6-1-20, Supplementary Data 2) to test whether a subrecent speleothem from this cave yields values comparable to modern temperatures. Samples were analyzed at 2.8 ± 0.2 ka and 4.4 ± 0.3 ka (Supplementary Data 2). Due to lack of water, we were unable to obtain data from younger layers of this stalagmite. The Holocene fluid inclusions yielded an average temperature of 2.0 ± 1.4°C, well within the error of the modern MAAT of 2.1°C. Further, temperatures in central Europe are reported to be within ±1°C of the mean temperature of the twentieth century at 2.5 ka 47 , suggesting that our subrecent fluid inclusion temperature provides a strong anchor point to LIG temperatures. Uncertainties reflect isotope measurement errors, δ 18 O/annual air temperature slope error, and one standard deviation of repeated measurements. The uncertainties are applied through all steps of the paleotemperature calculation. Further, uncertainties are propagated between SKR-FIT sampling locations.

Data availability
All data are available on Pangaea.