Speleothem Records from the Eastern Part of Europe and Turkey—Discussion on Stable Oxygen and Carbon Isotopes

: The region comprising of East Central Europe, South East Europe and Turkey contributed to the SISAL (Speleothem Isotopes Synthesis and AnaLysis) global database with stable carbon and oxygen isotope time-series from 18 speleothems from 14 caves. The currently available oldest record from the studied region is the ABA-2 ﬂowstone record (Abaliget Cave; Hungary) reaching back to MIS 6. The temporal distribution of the compiled 18 records from the region points out a ~20 kyr-long period, centering around 100 ka BP, lacking speleothem stable isotope data. The regional subset of SISAL_v1 records displays a continuous coverage for the past ~90 kyr for both δ 13 C and δ 18 O, with a mean temporal resolution of ~12 yr for the Holocene, and > 50 yr for the pre-Holocene period. The highest temporal resolution both for the Holocene and the pre-Holocene was achieved in the So-1 record (Sofular Cave; Turkey). The relationship between modern day precipitation δ 18 O (amount weighted annual and winter season mean values; 1961–2017) and climatological parameters was evaluated. The strong positive correlation found in East Central Europe reinforces the link between modern day precipitation δ 18 O, temperature and large-scale circulation (North Atlantic Oscillation) expected to be preserved in the speleothem δ 18 O record; while a negative relationship was documented between precipitation amount and oxygen isotope compositions in South East Europe. Variations of δ 13 C values are primarily interpreted as reﬂecting dry / wet periods across the region. Elevation gradients from three non-overlapping periods of the last ~5 kyr indicated elevation gradients around − 0.26% (cid:24) per 100 m − 1 for calcite δ 18 O.


Introduction
Southeastern Europe is the cradle of karst research [1][2][3]. The Kras Plateau in Slovenia gave its name to the entire discipline, and the very first karst research [4] and speleological [5] institutes were established in this region. Extended karstic terrains are found in various climatic zones, across the studied region, ranging from the extreme humid coastal areas of the Dinaric Karst (annual precipitation >4000 mm) to arid central Anatolia (annual precipitation <300 mm). Turkey and the southeastern part of Europe have a very diverse karst landscape, and a wide range of geochronological and geochemical investigations targeted the cave deposits in this region to infer information on past environments.  [50]. Numbers next to the sites refer to the site names in Table 1.
Only a few speleothem stable isotope records from the extended karst regions in the northern part of the studied are have been published (e.g., [17,56]. An explanation for this is that the permafrost conditions prevailed around [57] and beneath [58] the former Fennoscandian Ice Sheet, thereby preventing speleothem growth [17] because there was no water movement in the frozen ground. This hypothesis is indirectly supported by a recent study showing flowstone deposition in the Kraków-Częstochowa Upland at 975-470 ka BP [59]), predating the major glaciations in the region. The thawing of this frozen ground could be a prolonged process after deglaciation [60], maintaining unfavorable conditions for speleothem formation. The current karstwater table is near the surface in the Baltic region, creating phreatic conditions (i.e., completely submerged, or partly flooded) already in near surface cavities and subaerial speleothem formation has not been observed in known karst caves [61]. Speleothems were presumably formed in the Baltic region during the pre-Pleistocene [62], but the former vadose zone was denudated by extensive Pleistocene glaciation [63].
Additional identified records (e.g., from the Romanian Carpathians [64,65] the Low Tatras [66] or the Pontides [67][68][69]; Table 1) would help improving the spatiotemporal coverage of the region. A couple of records were not included to the SISAL database yet, either due to difficulties in obtaining data or due to the lack of essential metadata.  [50]. Numbers next to the sites refer to the site names in Table 1.
Only a few speleothem stable isotope records from the extended karst regions in the northern part of the studied are have been published (e.g., [17,56]. An explanation for this is that the permafrost conditions prevailed around [57] and beneath [58] the former Fennoscandian Ice Sheet, thereby preventing speleothem growth [17] because there was no water movement in the frozen ground. This hypothesis is indirectly supported by a recent study showing flowstone deposition in the Kraków-Częstochowa Upland at 975-470 ka BP [59]), predating the major glaciations in the region. The thawing of this frozen ground could be a prolonged process after deglaciation [60], maintaining unfavorable conditions for speleothem formation. The current karstwater table is near the surface in the Baltic region, creating phreatic conditions (i.e., completely submerged, or partly flooded) already in near surface cavities and subaerial speleothem formation has not been observed in known karst caves [61]. Speleothems were presumably formed in the Baltic region during the pre-Pleistocene [62], but the former vadose zone was denudated by extensive Pleistocene glaciation [63].
Additional identified records (e.g., from the Romanian Carpathians [64,65] the Low Tatras [66] or the Pontides [67][68][69]; Table 1) would help improving the spatiotemporal coverage of the region. A couple of records were not included to the SISAL database yet, either due to difficulties in obtaining data or due to the lack of essential metadata. Table 1. Metadata of the speleothem records from the Eastern part of Europe and Turkey. Column headings in italics represent field names that can be queried in the SISAL database. The numbers in parentheses next to the site names identify the site locations in Figure 1. Min/Max year correspond to the last/first ages of stable isotope data, respectively. BP is defined as years before present, where present is 1950 CE. Minimum and maximum age of stable isotope records were rounded to integer. Speleothems without an entity_id are not in SISAL_v1 [29,30].

Temporal Distribution of Speleothem Records
The oldest record in the region is the ABA_2 speleothem, collected from Abaliget Cave ( Figure 2; Table 1), which spans the period 160.6 to 140.2 ka BP. The period from~161 to~109 ka BP is represented by only four records with two overlapping segments. The first overlap occurs between the two Abaliget records from~143.5 to~140 ka BP, whereas the second overlap is seen between the Baradla records and the ABA_1 speleothem from 129 to~123 ka BP. The two Baradla records stopped growing at~109 ka BP, which is then followed by a ca. 20 kyr-long gap up to~89.7 ka BP ( Figure 2). From~90 ka BP onwards, there is continuous coverage of both δ 18 O and δ 13 C in the region, although the available number of isotopic records for the period 90-10 ka BP is low. The Dim-E3 record spans from~89.7 to~12.6 ka BP overlapping with the So-1 record after~50.3 ka BP (Table 1). These two speleothems are the only records from the region until 10.5 ka BP, when the Leany record starts ( Figure 2). The Holocene is the best represented period in the region (Figures 2-4), with the interval from~2.05 to~1.42 ka BP being covered by 7 records (PU-2, GK-09-02, MAR_L, POM2, CC-1, So-1 and S1). The speleothems Urşilor-and, Mavri Trypa caves' did not grow between~3 and~2.2 ka BP. Their restart in growth falls close to the start of the deposition of CC-1 from cave Ceremosnja at~2.4 ka BP. The previous observation of no-growth before cca. 2.3 ka in multiple speleothems ( Figure 4) might indicate a dry period in South East Europe, although the speleothems from caves Skala Marion, Kapsia and Ascunsă indicate otherwise.
The lack of data around 100 ka BP (from~109 ka BP to~89.7 ka BP in Figure 2) has been pointed out in a previous study [75], although radiometric ages from the Carpathian karstic region have been reported for a couple of speleothems during this period (e.g., [76,77]). Moreover, a recent compilation of radiometric ages from more than 90 speleothems and flowstones from the Croatian karst show an even temporal distribution back to the onset of MIS 5 [7], suggesting a good opportunity to develop speleothem stable isotope records for this period. In addition, there is a good chance to further extend the regional dataset beyond MIS 6, as records predating the ones in SISAL_v1 have already been identified from Turkish [69], Macedonian [51], Croatian [78] and Hungarian [79] localities. Moreover, calcareous cave deposits exceeding the current limit of applicability of U-Th dating method (~600 kyr) and for which U-Pb dating methods [80,81] could be applied, have been found in the region [59]. The average length of the records is ~13.8 kyr and the median is 4.9 kyr, including hiatuses. The shortest growth period was 54 yrs (POS-STM-4) and the longest ~77 kyr (Dim-E3). The temporal resolution of the SISAL stable isotope records ranges from sub-decadal (So-1, [38], S1 [42]; Trió, [11] to >1000 yr (Dim-E3, [44]). The average temporal resolution of the records is ~12 yr for the Holocene, while it is generally >50 yr for the pre-Holocene period ( Figure 4). The studies following the pioneering record of PU-2 in 2002 [47] achieved characteristically finer temporal resolution ( Figure  4), illustrating that both the sampling [84] and dating methods have highly improved over the past decades [85]. The average length of the records is~13.8 kyr and the median is 4.9 kyr, including hiatuses. The shortest growth period was 54 yrs (POS-STM-4) and the longest~77 kyr (Dim-E3). The temporal resolution of the SISAL stable isotope records ranges from sub-decadal (So-1, [38], S1 [42]; Trió, [11] to >1000 yr (Dim-E3, [44]). The average temporal resolution of the records is~12 yr for the Holocene, while it is generally >50 yr for the pre-Holocene period ( Figure 4). The studies following the pioneering record of PU-2 in 2002 [47] achieved characteristically finer temporal resolution (Figure 4), illustrating that both the sampling [84] and dating methods have highly improved over the past decades [85].

Availability of Cave Monitoring Data
Monitoring of climatic variables and stable isotope compositions of dripwater has been performed for at least one season in seven caves out of 14. In comparison to the whole SISAL_v1 database [30], the proportion of the monitored caves in the region is higher,~38% vs.~50%, respectively.
Monitoring of key in-cave variables (e.g., cave temperature, ventilation regime) influencing dripwater hydrology and carbonate saturation conditions provides the basis for a process-based understanding of modern environmental controls on stable isotope composition of speleothem and their growth conditions (e.g., [86][87][88]. In the lack of site-specific information (e.g., infiltration processes, calcite precipitation conditions), most of the paleoclimatic inferences of the early studies, e.g., [17,45], were not able to provide more than simplistic explanations. However, recently this knowledge-gap began to be filled as monitoring studies became abundant in the region [10,86,[89][90][91][92][93][94]. This predicts that much more experimental experience is being made available to help the interpretation of the speleothem recorded geochemical signatures over the region in the forthcoming studies.
For example, detailed monitoring enabled a recent study comparing modeled and measured dripwater δ 18 O values in Postojna Cave [95]. Monthly precipitation amount, δ 18 O data and monthly evapotranspiration were the input parameters, while mixing and delay processes taking place during infiltration were statistically simulated in the model. Authors found an excellent agreement between the measured δ 18 O spel data and the ones calculated using the modeled dripwater δ 18 O between 1984 and 2002. This argues that the effects and processes governing the δ 18 O spel variations are well understood, thus the δ 18 O spel signal in that particular system can reliably be interpreted.

Environmental Controls on Speleothem δ 18 O
The δ 18 O values in the surface waters of the Mediterranean Sea are higher than in the Black Sea, and much higher than in the Atlantic Ocean [96]. Together with the distinct evaporative conditions in the Mediterranean, this feature imprints the resultant atmospheric moisture with a source-specific isotopic signature, e.g., [97]. A key driver of changes in moisture transport trajectory in Europe during boreal winter (where most precipitation occurs) is the sea-level pressure (SLP) difference between the Arctic and the subtropical North Atlantic, the North Atlantic Oscillation (NAO; [98]). When the SLP difference is large (positive NAO mode), the main storm-track shifts to the north, transporting moisture and heat to western and northern Europe. When the SLP difference is weak (negative NAO mode), the North Atlantic jet is weakened and shifted to the south, bringing less heat, allowing an excess of Siberian cold air intrusions to the European continent, and bringing more moisture to southern Europe and its surroundings [99], including practically the entire study area. Besides the primary marine sources, recycled moisture from the Eurasian landmass at the eastern part of the region also becomes significant driver in continental areas [37,100,101]. Further, given the moisture transport trajectories from the ocean, weather systems travel a greater distance over land and are therefore affected by the processes of rainout, re-evaporation and moisture mixing in a positive NAO-mode. If precipitation is delivered from a more-distant moisture source (see above), then a larger fraction of the original moisture content is removed via subsequent precipitation events resulting in the continuous and progressive depletion of heavy isotopes in the remaining vapor, and lower δ 18 O (and also δ 2 H) in the fallen precipitation [102]. The source specific isotopic signature, modified by the previously mentioned driving processes imprinted into meteoric waters, is then transferred to cave dripwaters, providing the opportunity to use the water stable isotopes of speleothem archives (δ 18 O spel , δ 2 H of speleothem hosted fluid inclusions) as tracers for the intermittent moisture sources in different paleoclimatic conditions, assuming that the speleothem growth rate and the temporal persistence of atmospheric patterns bringing moisture from a given source renders it possible.
A recent global analysis on the oxygen isotope composition of cave drip water showed a pronounced seasonal bias toward the precipitation δ 18 O (δ 18 O prec ) of recharge periods in European caves located in a relatively warm and dry climate [103]. Based on the experience gathered from case studies of the study region, the stable isotope compositions of dripwater reflect the amount-weighted annual mean δ 18 O (and δ 2 H) of local precipitation, or display a winter dominance, depending on the site-specific infiltration conditions [90,92]. The results of a large-scale simulation model of karstic groundwater recharge over Europe and the Mediterranean have also pointed out that winter half year is the main period of effective infiltration in the study region [40]. Hence the effective δ 18 O/T slope relevant for the infiltrated water may be expected to fall between the annual and cold seasons' modern day δ 18 O/T slopes.
To conduct an analysis of modern day δ 18 O/T slopes in the region, monthly δ 18 O prec values were collected from the Global Network of Isotopes in Precipitation (GNIP [104]) stations of the study region along with monthly gridded precipitation amount, surface air temperature (E-OBS [105]) from the nearest grid cell to the GNIP stations. The time series of the Hurrell PC-Based NAO Index [106] were also gathered for the analysis. Amount-weighted annual and winter (October to March) mean δ 18 O values were computed. Only years for which GNIP δ 18 O prec data is available for enough months to represent >80% of the total annual precipitation and that had at least six years of evaluable data were used. The GNIP data provide a non-continuous coverage of 1961-2017; with an average of 23 and 24 years available for the annual and winter seasons respectively. Arithmetic means of air temperature ( • C), total amount of precipitation (mm) and NAO index were also computed for the same intervals. However, for the NAO index only the winter seasonal average was computed and used, because the annual aggregation of the NAO index may mask other atmospheric patterns [107]. The linear regression of the amount weighted δ 18 O prec data and meteorological parameters defined their long-term slopes/gradients. Spearman correlation coefficients (r s ) between δ 18 O prec data and the NAO index was used to evaluate the strength of the relationships.
As a result, a regional pattern may be observed from the modern day δ 18 O/T slopes (so-called "temperature effect"), making it possible to distinguish between the warm Mediterranean coastal sites and the continental ones. A progressive increase in the δ 18 O/T slope was detected as moving inland from the Mediterranean to Central Europe (Figure 5a). The temperature effect for annual precipitation ranges from −0.11% per • C for Mediterranean sites, and up to 0.75% per • C for continental locations (Figure 5a). For the winter half of the year, the values range between −0.04 and 0.83 % per • C, with the steepest gradients again found at the inland sites ( Figure 5a). As for the δ 18 O/precipitation, the relationship is mostly negative, and just occasionally significant in the Aegean region (p < 0.1), and at Zagreb in winter (p < 0.05) (Figure 5b).
According to additional observations, at lower latitudes (<45 • N) the δ 18 O prec variability negatively correlates with the winter NAO index, though to a rather insignificant degree, except for Ankara (r s = −0.43 p = 0.06; Figure 5c). On the contrary, north of the Adriatic Sea (>45 • N), the rank correlation relationship between δ 18 O and the NAO index is quite strong, and its significance increases towards the central parts of the continent in winter (Figure 5c).
The observed modern day δ 18 O/T slope is well above the temperature dependence of water-calcite isotopic equilibrium fractionation (−0.24% per • C [108]) at the continental locations of the Eastern part of Europe, while a different -insignificant-relation is observed at the coastal sites (Figure 5a). This δ 18 O/T difference implies that the oxygen isotope composition of speleothems from continental locations from the Eastern part of Europe may well reflect temperature variations, provided that the δ 18 O prec signal is transferred to the speleothem carbonate. As a supportive example, an excellent agreement has been found between speleothem δ 18 O spel of a modern stalagmite from Postojna Cave (Slovenia) and the regional record of O isotope composition of precipitation between 1984 and 2002 [95].
On the other hand, the lack of significant δ 18 O/T relationship, together with the occasionally negative correlation with the amount of precipitation in the southern part of the study region (Figure 5a,b) is in line with speleothem-based paleoclimate studies from the same area, which have found that δ 18 O spel strongly reflects fluctuations in precipitation amount [13,42,43]. Moreover, the correspondence between the interannual variability of total amount of October-January precipitation and δ 18 O spel from the Akcacale Cave (not yet in the database; Table 1) was sufficiently strong to allow the first quantitative reconstruction of amounts of precipitation in late fall and early winter [67,109].
A previous study assessed the spatial variability in European winter δ 18 O prec -NAO relationship using GNIP data available at that time [110]. No significant relationship was found for most of the stations analyzed in Eastern part of Europe and Turkey, and this study suggested a weak positive correlation (r s ) across this region [110]. However, the current re-evaluation, using an updated GNIP dataset has pointed out a significant and stronger relationship (r s difference > 0.35) between winter δ 18 O prec and NAO (avg. n = 24) across East Central Europe. A more recent research also documented positive correlations (δ 18 O prec and winter NAO) for Central Europe (GNIP sites west of 17 • E and north of 45 • N) [111], which are in agreement with the ones obtained from GNIP stations located in the western part of our study area (station codes D, E, F, H in Figure 5c). This provides a new insight into winter δ 18 O prec and NAO relationship in Europe.
Besides this continental sub-domain with significant positive correlation, the insignificant δ 18 O prec -NAO relationship is more restricted to the southernmost coastal areas (Figure 5c) than previously supposed. This spatial pattern revealed by the present analysis implies that the core region of the positive δ 18 O prec -NAO correlation mostly located in Central Europe [110] also extends eastwards. The significant negative correlation of δ 18 O prec -NAO in winter season found in Central Anatolia (station code A in Figure 5c) is an additional identified feature compared to the previous results [110].
This positive correlation in and around the Pannonian Basin, however, contradicts expectations. On the basis of the region's moisture source dynamics, in the case of a negative NAO mode, more Mediterranean 18 O-enriched (i.e., higher δ 18 O) moisture (see Gat et al. [97]) is to be expected at the inland locations, travelling via a shorter transport route that results in less negative δ 18 O prec . On the other hand, a positive NAO mode would deliver a dominantly 18 O-depleted (i.e., lower δ 18 O) Atlantic moisture to the continental locations along a relatively longer route that results in more negative δ 18 O prec . Since the observed relationship is contrary to the expected, other factors must be considered. A positive NAO state would not only be associated with Atlantic moisture dominance, but also with increased heat transport to Central Europe, with an inverse response (reduced heat transport) in Anatolia [99].
It follows from these considerations that the observed positive δ 18 O prec -T and δ 18 O prec -NAO relationships at the inland sites, along with the combination of positive δ 18 O prec -T and negative δ 18 O prec -NAO correlations in central Anatolia, can be attributed to the complex interplay of temperature-related oxygen isotope fractionation along the moisture transport route, plus the effect of moisture source composition. At the inland sites, temperature-related fractionation dominates over the source effect and the NAO influence appears via the heat transport. At the Anatolian site the moisture source effect exceeds the temperature-related isotope shift, inducing a negative δ 18 O prec -NAO correlation.

Environmental Controls on Speleothem δ 13 C
Based on the review of the speleothem-based paleoclimate studies from the region, it could be deduced that variations of δ 13 C values are primarily interpreted as a proxy for paleohydrological changes in the studied karstic areas although the highlighted effective/dominant responsible mechanism might be different. Low δ 13 C values indicate enhanced soil biological activity, which occurs in the case of a rather positive soil water balance (i.e., subsurface water recharge exceeds evapotranspiration) (Abaliget Cave [49], Ascunsă Cave [46], Baradla Cave [74], Trio Cave [11], Skala Marion Cave [43], Mavri Trypa Cave [42], Kapsia Cave [13]) or simplified to corresponding low (high) δ 13 C values with dry (wet) periods (Urşilor Cave [47], Modrič Cave [41]). The phosphorus contents of the Trió speleothem further support that the δ 13 C values reflect changes in soil biological activity induced by changes in soil hydro-balance [11]. Unlike in the majority of the area, in the warm semi-arid environments of the Anatolian karst, the secular changes in δ 13 C may track changes in the dominance of C3 vs. C4 in the vegetation type (Dim Cave [44], Sofular Cave [38]).
The situation can be more complicated for well ventilated caves, where degassing-related kinetic fractionation causes an increase in δ 13 C values, which partly overprints the paleohydrological signals. Two examples in the region are the Dim E3 and Leány records, which indicated enriched carbon isotope compositions (δ 13 C > −5% ) compared to the others in the region (Figure 2b) likely due to a ventilation effect; see also Section 5. However, minima in speleothem δ 13 C curves even in these cases can be regarded as largely unaffected by such kinetic processes [49] and likely reflect wet periods (i.e., positive soil water balance).

Regional Patterns in Carbon-and Oxygen Isotope Records
Prior to assessing the regional patterns, the records most prone to kinetic effect must be identified to avoid using records that may not be clear indicators of an environmental change [112]. A weak or no correlation between δ 13 C and δ 18 O spel data is conventionally indicative of deposition of calcite in isotopic equilibrium with the dripwater [113], therefore, the lack of significant correlation between the C and O isotope compositions, as in the case of the ABA_1 (p = 0.044), PU-2 (p = 0.251) and Dim-E4 (p = 0.089) records ( Figure S2, Table S1), suggest a quasi-equilibrium carbonate precipitation. In cases where the δ 13 C and δ 18 O values are positively correlated with a slope exceeding a value of 1, strong kinetic isotope fractionation effect can be expected to have occurred during carbonate precipitation [114]. The highest positive slope (2.64) with a significant correlation (adjusted R 2 = 0.45; p < 6.5 × 10 −15 ) between δ 13 C and δ 18 O spel values was observed in the CC-1 record. However, this is because the samples were gathered along a radial transect [45], rather than along the growth axis of the speleothem so that the kinetic effect is exaggerated as sampling points are closer to the edge. The strongest correlation (adjusted R 2 = 0.64) although with a less steep slope of 0.35 (p < 0.003) is seen in the short dataset of Pos-STM-4 record ( Figure S2, Table S1). This kinetic effect was not an issue in the original study, because its primary scope was to detect the increase in 14 C activity due to nuclear tests and to model the dynamics of soil carbon pools producing soil CO 2 [18]. In addition, other exposed records to kinetic effect were the Dim-E2, E3 and S1 records with significant positive δ 13 C-δ 18 O spel slopes (>1.4; p < 1.7 × 10 −15 ), which have to be handled with caution in advanced (paleo)climatological studies.
Nowadays, however, a weak correlation between δ 13 C and δ 18 O spel values is less taken to indicate isotopic equilibrium conditions. The replication test which investigates a high degree of coherence between individual δ 18 O and/or δ 13 C profiles from different speleothems from the same cave over the common period is a more stringent test of suitability for paleoclimatic studies [115]. Moving from the cave level (speleothems from the same cave [115]) to a regional level (speleothems from a cave system or nearby caves), the SISAL_v1 database and future versions provides a good opportunity for regional replication tests. On the one hand, if the records share common patterns, it might confirm a regional environmental signal, as in the case of e.g., South American Monsoon System [116]. On the other hand, trends and/or anomalies in stable isotope composition of speleothems occurring exclusively in an individual record may be owed to site specific factors. For example, the pronounced 13 C enrichment prevailing in the older section of the CC-1 record and preceding the growth stops of S1 is hardly replicated in the coeval speleothems from the region (Figure 3b). This suggests that at least some of their sections carry a strong kinetic effect and should be also treated with caution. It is worth mentioning, that these records (e.g., CC-1, S1) were also among the ones, possibly most affected by kinetic fractionation according to the δ 18 O spel vs. δ 13 C correlation.  [44]) speleothems, and inverse sudden changes in a north Anatolian (Sofular Cave [38]) record. The negative shift in So-1 resonates to a decreased heat transport from the cooled Atlantic inland. The mirrored pattern of δ 18 O from the Dim speleothems is in agreement with the Levantine records [117]. Hence, these variations can similarly be linked to hydroclimatological changes, i.e., less precipitation.
Regarding the full interval covered by the region's records, the negative peaks in δ 13 C at~120 ka BP and~43 ka BP are most profoundly replicated by two pairs of records: BAR-II#L and BAR-II#B and So-1 and Dim-E3 respectively (Figure 2b). However, the one at~43 ka BP is not associated with a concomitant 18 O peak ( Figure 2).
Comparing late-MIS 6 and MIS 2 glacial patterns, δ 18 O behave systematically (low δ 18 O spel values in glacial and high δ 18 O spel values in interglacial times), while the average δ 13 C of the ABA-1 record (Hungary) in the late-MIS 6 is lower and less variable than the average carbon isotope composition during MIS 2 in Anatolian records (Figure 2). This difference may derive from locally influenced carbon isotope values (e.g., related to ventilation and variable host rock contribution), but it may also be attributed to climatic conditions. Supporting the latter, paleoglaciological reconstructions argue for significantly more extended ice cover in the Dinarides and Hellenides in the MIS 6 compared to MIS 2 [118,119] suggesting more humid conditions. Paleoglaciological and paleobiological evidence argue for a similar MIS 2 vs. MIS 6 contrast regarding hydroclimate in the western Mediterranean (e.g., [120]). This hydroclimatological contrast might be partially reflected in speleothem isotopic records. Following the last glacial maximum (~20 ka BP) all available records (except Dim-E2) indicate an increase in δ 18 O spel . The highest change is displayed by the So-1 (~5% between 12 ka BP and 8 ka BP). This gradual drift in δ 18 O spel was linked to the long-term change in Black Sea moisture contribution [38].
In the Holocene, most of the records' δ 13 C values scatter around -10% and do not show any common trend or shift (Figure 3b). Important Holocene climate swings in Europe (e.g., 8.2 and 4.2 ka BP events) have been reported in some speleothem records of European subregions [31,121], however such records are not included in the regional-Eastern Europe and Turkey-subset of the current version of the SISAL_v1 database [30]. Negative δ 13 C anomalies around 3.5 ka BP were recorded in the Sofularand the Trio cave records, which were attributed to regional climate-change events [14]. The climate change starting at~3.7 ka BP with shifting to wetter conditions (lower δ 13 C values), was interrupted by a positive δ 13 C peak at~3.6 ka BP (Figure 3c). This prominent δ 13 C shift has been interpreted as a sign of environmental deterioration induced by the Thera volcanic eruption [14], which was later challenged on the basis of a significant time lag between the deposition of volcanic material and the positive δ 13 C shift [15].
The strong negative shift in δ 13 C around 3.6-3.5 ka BP is associated with a negative δ 18 O spel peak in both Trio-and So-1 records, likely related to a cold and wet phase [11]. Notably, the MAR_L record [43] replicates this δ 13 C pattern, although with different pattern in its δ 18 O values. The low sampling resolution of other currently available records covering the 4 to 3 ka BP period (PU-2, POM2; Figure 4) makes these records unsuitable for contributing to this issue. In the meanwhile, the high-resolution d-excess [122] record of the cave ice deposit from Scărişoara Ice Cave (Romania) shows an abrupt increase at that time, likely suggesting more Mediterranean moisture contribution, accompanying wetter conditions [21].

Assessment of Spatial Patterns of Speleothem δ 18 O in the Holocene
The SISAL_v1 database contains 12 records from the Eastern part of Europe and Turkey covering the past 11.7 kyr, with δ 18 O spel values ranging between -3 and -12% (Figure 3a). The mean δ 18 O values tend to follow a latitudinal trend, albeit, two records stand out from this pattern. MOD-22 has the highest δ 18 O spel values despite not being at the southernmost location, while the So-1 has relatively low δ δ 18 O spel values (Figure 3a) despite its host cave not being the northernmost one (Figure 1; Table 1), thus contradicting a putative latitudinal gradient. These low values of So-1 could be driven by moisture, derived from the less evaporated nearby Black Sea [38,39], compared to the strongly evaporated surface waters of the Mediterranean Sea feeding other areas of the region [36].
The relatively higher abundance of data from the mid-Holocene (Figure 3a) provides an additional insight into the spatial trends of δ 18 O spel values. Average δ 18 O spel values were calculated for contemporarily formed deposits from three separate semi-millennial periods from the Holocene ( Figure 6). Unlike previous observations [123] continental trends could not be observed probably due to the limited spatial extent of the available data (latitudinal extent <12 • , longitudinal extent <18 • ; Figure 1). Additionally, the infiltration of precipitation water at individual cave sites is seasonally biased [90,92], resulting in drip waters with varying winter/summer precipitation contributions and consequently changing δ 18 O values.
Quaternary 2018, 2, x FOR PEER REVIEW 18 of 26 Figure 6. Elevation gradient of δ 18 Ospel for three semi-millennial time horizons for the past 5000 yrs for SE Europe. The italic numbers above the dots refer to the site names found in Table 1.
However, we found negative correlations between δ 18 Ospel records from eastern and southeastern Europe and site elevation ( Figure 6). This stands even though "elevation m amsl" (Table  1) in the database corresponds to the cave entrance elevation [29,30] rather than the elevation of dripwater recharge area, the latter being the actual controlling factor. It can be assumed that the difference in site elevation primarily represents the difference in infiltration area elevations.
Despite this uncertainty regarding the elevation of the actual recharge area, estimated δ 18 Ospel elevation gradients for SE Europe calculated for three non-overlapping semi-millennial periods (and excluding records with high risk of kinetic effects as pointed out in Section 5) during the past ~5 kyr suggests that this elevation trend is characteristic of the Holocene. Computed elevation gradients were significant (at p < 0.1 level) in all but one cases and scatter around −0.26‰/100 m.  Table 1. However, we found negative correlations between δ 18 O spel records from eastern and southeastern Europe and site elevation ( Figure 6). This stands even though "elevation m amsl" (Table 1) in the database corresponds to the cave entrance elevation [29,30] rather than the elevation of dripwater recharge area, the latter being the actual controlling factor. It can be assumed that the difference in site elevation primarily represents the difference in infiltration area elevations.
Despite this uncertainty regarding the elevation of the actual recharge area, estimated δ 18 O spel elevation gradients for SE Europe calculated for three non-overlapping semi-millennial periods (and excluding records with high risk of kinetic effects as pointed out in Section 5) during the past~5 kyr suggests that this elevation trend is characteristic of the Holocene. Computed elevation gradients were significant (at p < 0.1 level) in all but one cases and scatter around −0.26% /100 m.
The documented negative shift in δ 18 O spel values with increasing elevation is regulated by two main factors. Depletion of meteoric precipitation in heavy isotopes with elevation is expected [124] and isotopic "elevation effect" in precipitation has been measured in this region as ranging from −0.24% per 100 m to −0.39% per 100 m [37,88,[125][126][127] transferring with a less steep gradient to shallow groundwater (e.g., −0.24 % per 100 m, [128], interquartile margins: −0.11 and −0.24% per 100 m, ref. [129] and dripwater (e.g., −0.11% per 100 m, [88]). This process may account for most of the variance. Mean surface temperature also decreases with elevation, which is reflected in cave air temperature. Temperature dependent calcite-water fractionation-with a gradient of −0.24% per • C [108]-likely plays an additional role.
Whether the δ 18 O spel -elevation gradient ( Figure 6) prevails over geological timescales, especially under distinct climate settings such as the LGM, can be better explored with an increased number of records, providing a much-improved replication in this specific period. The robustness of these empirical gradients can be further explored extending the study region towards nearby regions [31,117].

Stable Isotopic Composition of Fluid Inclusions of Speleothems from South East Europe and Turkey to Complement Paleoenvironmental Implications of δ 18 O spel
Stable isotope analyses of inclusion-hosted water were performed on about a third of the speleothems records (n = 18) included in the SISAL_v1 database from the region; six speleothems belonging to four caves. Note, that the fluid inclusion data of the ABA 1 and ABA 2 flowstone cores [49] are unpublished. Out of the 371 speleothem records in the entire SISAL_v1 database, the isotopic composition of fluid inclusions was analyzed only in 27; 6 out of 18 in the study region. Despite the relatively small number of locations studied, stable isotope composition of inclusion-hosted water as a complementary parameter is well represented in the region.
This regional database provides an opportunity to illustrate the auxiliary information of fluid inclusion δ 2 H in interpreting carbonate stable isotope compositions for paleoclimate inferences. The sub-regional Middle Bronze Age (~3.5 ka BP) climate swing recorded by certain Hungarian, Turkish and Greek speleothems' δ 13 C and δ 18 O records (Section 5.1) coincides with a~15% negative shift of fluid inclusion δ 2 H of the Trió speleothem peaking at 3.5 ka BP (Figure 3a). This was interpreted as a combined effect of temperature change and an increased contribution of winter precipitation to the infiltration [11]. Detailed paleoclimatological reconstruction for this period is important because archaeological evidence points to substantial societal-changes in the European-Mediterranean region occurring at the same time [11,130,131]. A similar series of temperature and precipitation seasonality changes were inferred for the last interglacial on the base of coupled negative shifts in δ 13 C, δ 18 O spel and fluid inclusion δ 2 H [48]. These complex paleoenvironmental changes could not have been revealed without the combined assessment of carbonate and fluid inclusion stable isotope compositions illustrating that latter serves as an important paleotemperature proxy.

Conclusions and Outlook
The diverse karst landscapes of the Eastern part of Europe and Turkey offer a great potential for reconstructing past environmental changes based on isotope records from speleothems. The regional subset of the first global database on carbonate stable isotope composition (SISAL_v1 [29,30]) includes 18 speleothem records. The number of overlapping records in this subset is generally low, with the most represented period spanning from 5-1 ka BP with 5 to 8 records. The Sofular Cave So-1 record is the only one in the region covering the last glacial cycle and extending into the Holocene. This is an obvious weakness in temporal coverage that should to be improved in the future. Thus, there is a strong need to gather speleothems from this region that provide a continuous record of the environmental response induced by a transition from cold glacial to warm interglacial climate regime.
As documented by most of the case-studies, stable carbon isotope compositions of speleothems of the region generally reflect hydroclimatological signals via related changes in soil activity. In the present study, a spatially complex picture was observed when the modern-day δ 18 O prec data were compared with annual and winter climate variables (temperature, precipitation amount and NAO) across the Eastern part of Europe and Turkey ( Figure 5). The generally positive δ 18 O-T relationships are associated with positive δ 18 O-NAO correlations in eastern Central Europe, and a statistically significant negative δ 18 O-NAO coefficient in Central Anatolia. These relationships can be attributed to the complex interplay of temperature-dependent oxygen isotope fractionation in meteoric waters and δ 18 O differences in the moisture source regions, i.e., 18 O-depleted Atlantic and 18 O-enriched Mediterranean areas. The strong positive δ 18 O-T relationships at the inland sites reflect the effect of NAO-driven heat transport that should be evaluated at the individual locations of future paleoclimatic studies.
The limited number of records in the regional subset of SISAL_v1 reflect the expected contrast between glacial and inter-glacial conditions. Low δ 18 O spel values dominated glacial periods and high ones characterized interglacial times. Regarding glacial times, the few available glacial δ 13 C records (Abaliget Cave (S Hungary), Sofular-and Dim caves (Turkey)) suggest a difference between MIS6 and MIS2 periods with the former assumed to be more humid than the latter.
Systematic δ 13 C or δ 18 O spel shifts or drifts have not been observed in the Holocene apart from the period of 4 to 3 ka BP, which is characterized by a series of corresponding isotope shifts in different records. This presumed Middle Bronze Age climate fluctuation is further confirmed by speleothem derived fluid inclusion δ 2 H from the area, which aid the complex paleoclimatological interpretation [11]. It would be intriguing to further explore the presence of the "~3.5 ka event" in other records from the Eastern part of Europe and Turkey and compare the observations on a larger geographical scale.
The limited horizontal extent of the area covered by the regional subset, did not allow to confirm the existence of longitudinal trends in δ 18 O spel data across south-eastern Europe during the Holocene. However, examination of elevation gradients from South East Europe suggested a long-term elevational gradient in δ 18 O spel values (around −0.26% 100 m −1 ) during the Holocene.
Expansion of the database by including recently published data, and exploration and exploitation of regions such as the Dinaric Karst [7] will yield more information in the near future. The number of identified, but not yet included records (Table 1) will obviously increase the probability to explore common temporal patterns in δ 13 C and δ 18 O spel data in a spatially robust way. The combined use of additional parameters complementary to carbonate stable isotopes from speleothems (e.g., radiocarbon activity data, trace element content, δ 2 H of fluid inclusion) and other sources (e.g., paleobiological and archaeological information, etc.) are needed to improve our understanding on the environmental response of surface-and subsurface systems to climate changes. Lastly, it should be emphasized that speleothem derived paleoenvironmental data-typically unevenly spaced and associated with timescale error should be assessed using proper statistical tools [132,133] capable to provide more robust results in speleothem paleoclimatology than before.
Supplementary Materials: The following are available online at http://www.mdpi.com/2571-550X/2/3/31/s1; Figure S1: Climate of the Eastern Europe & Turkey SISAL region; Figure S2: Cross plot of δ 13 C vs. δ 18 O; Table S1: Statistics for the regression models δ 13 C vs. δ 18 O. Acknowledgments: Authors thank all colleagues involved in SISAL for intriguing discussions and collaboration on the preparation of this manuscript. SISAL is a working group of the Past Global Changes (PAGES) programme and we thank PAGES for their support of this activity. We are grateful to Andy Baker, Sandy Harrison and Laia Comas-Bru for detailed and constructive criticism. Authors also thank the members of the World Karst Aquifer Mapping project for providing the shape file used to derive the carbonate aquifer maps ( Figure 2) and acknowledge the E-OBS dataset from the EU-FP6 project UERRA (http://www.uerra.eu) and the Copernicus Climate Change Service, and the data providers in the ECA&D project (https://www.ecad.eu). In addition, we thank Laia Comas-Bru for plotting Figures 1 and 4 for the study.