Reconstructing seasonality through stable-isotope and trace-element analyses of the Proserpine stalagmite, Han-sur-Lesse cave, Belgium: indications for climate-driven changes during the last 400 years

Fast-growing speleothems allow for the reconstruction of palaeoclimate down to a seasonal scale. Additionally, annual lamination in some of these speleothems yields highly accurate age models for these palaeoclimate records, making these speleothems valuable archives for terrestrial climate. In this study, an annually laminated stalagmite from the Han-sur-Lesse cave (Belgium) is used to study the expression of the seasonal cycle in northwestern Europe during the Little Ice Age. More specifically, two historical 12-year-long growth periods (ca. 1593–1605 CE and 1635–1646 CE) and one modern growth period (1960–2010 CE) are analysed on a sub-annual scale for their stable-isotope ratios (δ¹³C and δ¹⁸O) and trace-element (Mg, Sr, Ba, Zn, Y, Pb, U) contents. Seasonal variability in these proxies is confirmed with frequency analysis. Zn, Y and Pb show distinct annual peaks in all three investigated periods related to annual flushing of the soil during winter. A strong seasonal in-phase relationship between Mg, Sr and Ba in the modern growth period reflects a substantial influence of enhanced prior calcite precipitation (PCP). In particular, PCP occurs during summers when recharge of the epikarst is low. This is also evidenced by earlier observations of increased δ¹³C values during summer. In the 17th century intervals, there is a distinct antiphase relationship between Mg, Sr and Ba, suggesting that processes other than PCP, i.e. varying degrees of incongruent dissolution of dolomite, eventually related to changes in soil activity and/or land-use change are more dominant. The processes controlling seasonal variations in Mg, Sr and Ba in the speleothem appear to change between the 17th century and 1960–2010 CE. The Zn, Y, Pb, and U concentration profiles; stable-isotope ratios; and morphology of the speleothem laminae all point towards increased seasonal amplitude in cave hydrology. Higher seasonal peaks in soil-derived elements (e.g. Zn and Y) and lower concentrations of host-rock-derived elements (e.g. Mg, Sr, Ba) point towards lower residence times in the epikarst and higher flushing rates during the 17th century. These observations reflect an increase in water excess above the cave and recharge of the epikarst, due to a combination of lower summer temperatures and increased winter precipitation during the 17th century. This study indicates that the transfer function controlling Mg, Sr and Ba seasonal variability varies over time. Which process is dominant – either PCP, soil activity or dolomite dissolution – is clearly climate driven and can itself be used as a palaeoenvironment proxy.

1 Introduction

Speleothems have been successfully used to reconstruct palaeoclimate on various timescales (Fairchild and Baker, 2012), from tropical latitudes (e.g. Wang et al., 2001) to temperate areas (e.g. Genty et al., 2003). Fast-growing speleothems enable palaeoclimate reconstructions to reach seasonal resolution, or even higher (Van Rampelbergh et al., 2014). Some speleothems are characterized by distinct physical and/or geochemical layering, which improves chronologies and lends more confidence to the interpretation of proxy records at these high temporal resolutions (e.g. Mattey et al., 2008; Boch and Spötl, 2008). The occurrence of visible annual laminae in speleothems has been reported from sites all over the world (Baker et al., 2008). A common expression of this visible layering is an alternation of dark compact laminae (DCL) and white porous laminae (WPL), as defined by Genty and Quinif (1996). According to Dreybrodt (1999) and Baker et al. (2008), the origin of visible seasonal layering is related to seasonal variations in drip rate and in drip water supersaturation. Additionally, seasonal changes in cave ventilation can influence the intensity of the degassing process of the drip water and influence the pH of the fluid from which the speleothem precipitates, producing seasonal variations in crystal fabric (e.g. Frisia et al., 2003). These two processes can work in conjunction, with drip rate influencing annual laminae thickness while degassing influences the speleothem fabric (Boch et al., 2010). Such changes in drip rate often coincide with the presence of a varying degree of prior calcite precipitation (PCP). PCP is the process of calcite precipitation upstream of the site of speleothem deposition (Fairchild et al., 2000). An increase in PCP occurs when the ability of infiltration waters to degas increases. Therefore, a higher degree of PCP is attributed to drier periods (Fairchild et al., 2000; Fairchild and Treble, 2009). Variations in the amount of PCP have been observed on a seasonal scale (e.g. Johnson et al., 2006).

The presence of seasonally laminated speleothems in Belgian cave systems has been known for several decades (e.g. Genty and Quinif, 1996). The best known example is the Proserpine stalagmite, which was cored in the cave of Han-sur-Lesse and first studied by Verheyden et al. (2006). The speleothem has a well-expressed visual and geochemical seasonal layering over the last 500 years as inferred from and U ∕ Th dating and lamina counting between the radiometric ages (Van Rampelbergh et al., 2015). This geochemical layering is reflected by sub-annual variations in stable-isotope ratios (δ¹³C and δ¹⁸O). A thorough understanding of modern seasonal control on variations in δ¹³C and δ¹⁸O in speleothem calcium carbonate results from rigorous monitoring of the conditions at the sample site in Han-sur-Lesse cave as carried out by Van Rampelbergh et al. (2014) for the period 2012–2014.

In addition to the commonly used speleothem δ¹⁸O and δ¹³C proxies, the use of trace-element concentrations (e.g. Mg, Sr, Ba, Zn and U) as palaeoclimate and palaeoenvironmental proxies is becoming standard practice in speleothem reconstructions (Fairchild et al., 2000; Treble et al., 2003; Fairchild and Treble, 2009; Griffiths et al., 2010; Regattieri et al., 2016). The use of trace elements brings additional information that can be used to unravel seasonal variability in speleothem chemistry. Examples of this include the use of trace-element concentrations as proxies for precipitation (Wang et al., 2001; Baldini et al., 2002), soil processes (Regattieri et al., 2016) or changes in sediment supply (Regattieri et al., 2016). They have also been linked to changes in the atmospheric load of anthropogenic and volcanic-derived aerosols, as well as volcanic ash fall events from speleothem records (Frisia et al., 2005; Jamieson et al., 2015).

The first objective of this study is to better characterize the geochemical layering of the Proserpine speleothem by adding trace-element proxies to improve the understanding of processes driving the geochemical layering and to further resolve its relation with seasonal climatic variability. In addition, the present study compares the seasonal cycle within earlier identified cold periods (Verheyden et al., 2006; Van Rampelbergh et al., 2015; Fig. S1 in the Supplement) to present-day seasonal signals. To achieve this, two 12-year-long stalagmite growth periods (1593–1605 CE ± 30, hereafter P16, and 1635–1646 CE ± 30, hereafter P17) and a more recent growth period deposited between 1960–2010 CE (hereafter referred to as P20) were analysed at a sub-annual scale for their stable-isotope (δ¹³C and δ¹⁸O) and trace-element variations. The names of these growth periods (roughly) refer to the century in which they were deposited (16th, 17th and 20th century, respectively). Chemical changes in these three intervals will be interpreted in the context of the longer stable-isotope record that is available from Van Rampelbergh et al. (2015) to discuss changes in hydroclimate and seasonality over the past 500 years as recorded in the Proserpine speleothem.

Figure 1 (a) Location of the Han-sur-Lesse cave system (50^∘06^′51^′′ N, 5^∘12^′12^′′ E) with the entrance and exit of the Lesse River, the Salle-du-Dôme and the Père-Noël Cave. North is upwards. (b) Map showing the location of the Proserpine stalagmite within the Salle-du-Dôme. The insert shows the position of the core retrieved from the speleothem. Images adapted from Van Rampelbergh et al. (2014, 2015). (c) Lithological column of the Devonian strata in the Synclinorium of Dinant including the two formations (Fromellennes and Mont d'Haurs) in which the Han-sur-Lesse cave is situated. Note the presence of partly dolomitized biostromal limestones and shales in both these formations.

2 Geological setting

2.1 Han-sur-Lesse cave

With a total length of approximately 10 km, the Han-sur-Lesse cave system, located within a limestone belt of Middle Devonian age, is the largest known subterranean karst network in Belgium (Fig. 1a). The cave system was formed by a meander cut-off of the Lesse River within the Massif de Boine, which is part of an anticline structure consisting of Middle-to-Late Givetian reefal limestones (i.e. the Mont-d'Haurs and Fromelennes formations (Fm); Delvaux de Fenffe, 1985; Fig. 1c). The limestone epikarst reaches a thickness of 20–50 m above the cave, as estimated by the map of the cave and the surface. (Quinif, 1988). Studies have shown the presence of dolomite in these Givetian limestones which are directly overlying the Han-sur-Lesse cave (Fig. 1). Within the Mont-d'Haurs Fm., the biostromal limestones are alternated with fine-grained micritic limestones and dolomitic shales (Preat et al., 2006). Additionally, studies by Verheyden et al. (2000) and Pas et al. (2016) on the Middle Devonian strata in which the Han-sur-Lesse cave is situated have shown that dolomitized beds also occur within the limestones of the Fromelennes Fm. Specifically, the Salle-du-Dôme, in which the Proserpine is located, opens in the Devonian Givetian limestone beds forming an anticline structure, which explains the surface geomorphology of the hill in which the cave is located (Delvaux de Fenffe, 1985; Fig. 1). Since no impermeable formations are present above the cave, precipitation directly seeps through the thin (∼25 cm) soil and enters the epikarst.

Figure 2 Chart showing the calculated theoretical amount of water excess calculated with the Thornthwaite equation (Thornthwaite and Mather, 1957), based on temperature and precipitation data near Han-sur-Lesse cave from 1999 to 2012 (Royal Meteorological Institute, KMI). X axis represents the months from January (1) to December (12).

The Han-sur-Lesse cave is located ∼200 km inland at an elevation of 200 m a.s.l. The region is marked by a warm temperate, fully humid climate with cool summers, following the Köppen–Geiger classification (Kottek et al., 2006). In the period 1999–2012, the region experienced a seasonality in monthly temperatures between 2.5 ^∘C (January) and 18 ^∘C (July), with a mean annual temperature of 10.2 ^∘C (Fig. 2). Precipitation is fairly constant year-round around 40–50 mm month⁻¹, with 2 months (July and August) experiencing increased precipitation of 90–95 mm month⁻¹. Mean annual precipitation was 820 mm yr⁻¹ in Rochefort, 10 km from Han-sur-Lesse (Royal Belgian Meteorological Institute, Brussels, Belgium). The study site is affected by a North Atlantic moisture source all year round (Gimeno et al., 2010) and the amount of precipitation does not follow a seasonal distribution. Calculations applying the Thornthwaite formula (Thornthwaite and Mater, 1957) show that there is a strong seasonal trend in the water excess, i.e. the amount of rainfall minus the amount lost by evapotranspiration, with water excess only occurring from October to April (Genty and Deflandre, 1998; Genty and Quinif, 1996).

The studied speleothem was retrieved from the Salle-du-Dôme in 2001, a 150 m wide and 60 m high chamber that formed by collapse of the roof of the limestone (Fig. 1b). The Salle-du-Dôme is well ventilated, as it is located close to the cave exit and connected through two passages to nearby chambers. Monitoring of cave atmosphere within the Salle-du-Dôme for the period 2012–2014 showed that in 2013 the temperature inside the chamber varied seasonally between 10.5 and 14.5 ^∘C (Van Rampelbergh et al., 2014). Similar seasonal trends in temperature are observed for drip water sampled at the site of the Proserpine speleothem, but the mean is 0.5 ^∘C colder than the outside air temperature. The pCO₂ of the cave air fluctuates between 400 and 1000 ppmv (in July and August), and it averages 500 ppmv during the year. In summer (July–August), a rapid and temporary (2-month) increase to 1000 ppmv is observed. Also during summer, rainwater δ¹⁸O and δD above the cave increase by 3 ‰ and 30 ‰ (VSMOW, Vienna Standard Mean Ocean Water) respectively, likely due to the atmospheric temperature effect as described by Rozanski et al. (1992). In contrast, drip water δ¹⁸O and δD remain fairly stable throughout the year, with means of −7.65 ‰ and −50.1 ‰ VSMOW and standard deviations of 0.07 ‰ and 0.6 ‰ VSMOW, respectively. During late summer (September), an increase of 1.5 ‰ is observed in the δ¹³C values of dissolved inorganic carbon (DIC) within the drip water

2.2 Proserpine speleothem

The Proserpine speleothem is a 2 m high, stalagmitic flowstone with a flat top. The speleothem has a surface area of 1.77 m² and is fed by a drip flow with drip rates ranging between 100 and 300 mL min⁻¹. The type of speleothem is described in the literature as a “tam-tam speleothem” and is characterized by the occurrence of millimetre-to-centimetre-scale rimmed pools, which causes chaotic “pool-like” deposits to occur in parts of the speleothem (Verheyden et al., 2006). The speleothem grew over a period of approximately 2 kyr and thus has an exceptionally high mean growth rate of 1 mm yr⁻¹. This large speleothem was drilled and a 2 m long core was retrieved. The upper 50 cm of this core, dating back to approximately 1500 CE (Fig. S2), shows a well-expressed layering of alternating dark compact layers (DCL) and white porous layers (WPL; Verheyden et al., 2006). Previous studies concluded that simultaneous multi-decadal changes in different proxies (such as crystal fabric, growth rate, layer thickness, and oxygen and carbon stable-isotope ratios) indicate that these are controlled by climatic, environmental or anthropogenic factors and that some parts of the Proserpine speleothem have been deposited out of isotopic equilibrium with the drip water (Verheyden et al., 2006; Van Rampelbergh et al., 2015).

Based on a detailed cave-monitoring study at the Proserpine site in the years 2012 to 2014, Van Rampelbergh et al. (2014) showed that δ¹⁸O and δ¹³C of seasonally deposited calcite reflect isotopic equilibrium conditions and that variations in stable-isotope ratios are induced by seasonal changes. These seasonal changes in stable-isotope ratios correspond with the observed visible layering. The speleothem δ¹⁸O value is believed to reflect changes in seasonal cave climatology. While drip water δ¹⁸O remains constant, calcite δ¹⁸O decreases by ∼0.6 ‰ in summer months, caused by temperature-dependent fractionation during calcite precipitation. This fractionation was calculated to be −0.2 ‰ ^∘C⁻¹ (Van Rampelbergh et al., 2014). In contrast, δ¹³C reflects seasonal changes occurring at the epikarst level. A ∼1.5 ‰ increase in δ¹³C in drip water DIC during late summer is directly reflected in the freshly deposited speleothem calcite. The enrichment in drip water δ¹³C values occurs shortly after the observed decrease in drip water discharge, and therefore seasonal variations in the degree of prior calcite precipitation in the epikarst has been hypothesized to be the main driver of seasonal δ¹³C changes in the drip water (Van Rampelbergh et al., 2014).

Table 1 Schematic overview providing the observed changes and interpretation for the different proxies of P20, P17 and P16. PCP represents prior calcite precipitation, IDD represents incongruent dissolution of dolomite, DCL represents dark compact layers, and WPL represents white porous layers.

2.3 Dating

The age–depth model of the Proserpine speleothem core has been established and discussed by Van Rampelbergh et al. (2015) and is provided in the Supplement (Fig. S2). This age–depth model was constructed by using a combined approach of U-Th radiometric dating, based on 20 U-Th ages, and layer counting. It was shown that the number of counted layers is in good agreement with the U-Th ages (see Table 1 in Van Rampelbergh et al., 2015). However, 9 to 10 cm from the top of the core, a perturbation with heavily disturbed calcite occurs, making it impossible to construct a continuous layer-counting chronology. Remains of straw and soot were found within this perturbation, suggesting that at that time, fires were lit on the speleothem's palaeosurface (Verheyden et al., 2006). Layer counting gave an age of 1857 ± 6 CE for the reestablishment of calcite deposition after the perturbation and U-Th age–depth modelling showed that the start of the perturbation occurred at 1810 ± 45 CE (Van Rampelbergh et al., 2015). Radiocarbon dating of the straw fragments embedded in the calcite gave an age between 1760 and 1810 CE, with 95.4 % probability. The age of 1810 ± 45 CE is used to restart the layer counting after the perturbation towards the bottom of the core. This gave an age of 1593 to 1605 ± 30 CE for the P16 and 1635 to 1646 ± 30 CE for the P17 section. The more recent section P20 studied here is situated above the perturbation and its age could be confidently established through annual layer counting between 1960 and 2001 CE. The same chronology for the Proserpine speleothem was previously used in conjunction with petrographic evidence to show that decreases in cave visitation coincided with the two World Wars, highlighting the accuracy of the age model (Verheyden et al., 2006).

Figure 3 (a) Overview of a long (∼500 yr) record of stable-isotope ratios and annual layer thickness through the Proserpine speleothem measured by Van Rampelbergh et al. (2015). Red boxes indicate the locations of high-resolution transects discussed in this study. (b) The three studied growth periods P20 (1960–2010 CE), P17 (1633–1644 ± 30 CE) and P16 (1593–1605 ± 30 CE) are shown. The yellow rectangles mark the sections that were drilled and sampled for δ¹³C and δ¹⁸O analysis; the red lines represent the LA-ICP-MS transects. Numbers in grey indicate the observed layer couplets. Note that the images (b) of parts of the Proserpine speleothem that contain the three growth periods are oriented vertically (top of the picture upwards) but that the transects themselves are oriented parallel to the local growth direction (perpendicular to the growth laminae marked in blue). This results in transect P17's orientation to be at a ∼25^∘ angle with respect to the image (and the vertical).

3 Methods

3.1 Analytical procedures

The three growth periods studied are shown in Fig. 3 and their age is derived from an age–depth model based on U-Th dating and layer counting (Verheyden et al., 2006; Van Rampelbergh et al., 2015; Fig. S5; see Sect. 2.3). The number of counted layer couplets over the annually layered 500 years determined the seasonal character of the layers and demonstrated that two layers (one DCL and a WPL) were deposited per year. The number of years obtained by layer counting between two U-Th datings was compared with the number of years suggested by the U-Th ages. We combine results of both independent dating methods to produce the final age model. All growth transects are sampled parallel to the local direction of growth of the Proserpine speleothem and perpendicular to the growth laminae (see Fig. 2).

For δ¹³C and δ¹⁸O analysis, powder samples were acquired with a Merchantek micromill (Merchantek, Electro Scientific Industries Inc. (ESI), Portland, OR, USA, coupled to a Leica GZ6, Leica Microsystems GmbH, Wetzlar, Germany) equipped with tungsten carbide dental drills with a drill bit diameter of 300 µm. The powders are stored in a 50 ^∘C oven prior to the analysis to avoid δ¹³C and δ¹⁸O isotopic contamination. Measurements for P16 and P17 are carried out on a Nu Perspective isotope ratio mass spectrometer (IRMS) coupled to a Nucarb automated carbonate preparation device (Nu Instruments, UK) at the Vrije Universiteit Brussel (Belgium). The δ¹³C and δ¹⁸O records of P16 and P17 consist of 201 and 116 data points, respectively, resulting in temporal resolutions of ∼20 and ∼10 data points per year, respectively. The analysis of the P20 interval is an extension of the previously published seasonally resolved 1976–1985 transect (Van Rampelbergh et al., 2014) and was carried out on a Delta Plus XL IRMS coupled to a Kiel III carbonate preparation unit (Thermo Fisher Scientific, Germany) also at the Vrije Universiteit Brussel. For P20, a total of 350 samples were analysed, providing a temporal resolution of approximately seven data points per year. All results are displayed as ‰ VPDB (Vienna Pee Dee Belemnite) with the individual reproducibility reported as 2 standard deviation (SD) uncertainties. Within each batch of 10 samples, the in-house reference material MAR2-2, prepared from Marbella limestone and calibrated against NBS-19 (Friedman et al., 1982), is measured together with the samples to correct for instrumental drift (δ¹³C: 3.41±0.10 ‰ (2 SD) VPDB; δ¹⁸O: 0.13±0.20 ‰ (2 SD) VPDB). Total uncertainties for δ¹³C and δ¹⁸O are 0.03 ‰ and 0.09 ‰ (1 SD), respectively, for the Nu Perspective setup. With the Delta Plus XL setup these are slightly higher, being 0.04 ‰ and 0.10 ‰ for δ¹³C and δ¹⁸O, respectively (Van Rampelbergh et al., 2014).

Trace-element variations are determined using inductively coupled plasma mass spectrometry complemented by a laser ablation sample introduction system (LA-ICP-MS) at Ghent University (Belgium). The LA-ICP-MS setup consists of a 193 nm ArF excimer Analyte G2 laser ablation system (Teledyne Photon Machines, Bozeman, MT, USA) coupled to a single-collector sector field “Element XR” ICP-MS unit (Thermo Fisher Scientific, Bremen, Germany). The laser was used to sample adjacent positions along a line segment parallel to the growth axis. The positions were ablated one by one for 15 s with a laser spot size of 50 µm in diameter, a repetition rate of 30 Hz and a beam energy density of 3.51 J cm⁻². The line segments for P16, P17 and P20 are drilled at 287, 249 and 445 individual positions, respectively. Sampling via individual laser spots was preferred over the conventional approach of continuous line scanning because the single positions can be sampled longer, resulting in an improved limit of detection. To carry out the analyses, the speleothem sections and reference materials were mounted in a HELEX 2 double-volume ablation cell. The helium carrier gas (0.5 L min⁻¹) is mixed with argon make-up gas (0.9 L min⁻¹) downstream of the ablation cell, and introduced into the ICP-MS unit, operated in low mass-resolution mode. Transient signals for magnesium (Mg), aluminium (Al), silicon (Si), phosphorus (P), sulfur (S), potassium (K), iron (Fe), manganese (Mn), zinc (Zn), rubidium (Rb), strontium (Sr), yttrium (Y), barium (Ba), lead (Pb), thorium (Th) and uranium (U) were monitored during analysis of the laser-induced aerosol. Cool plasma conditions (800 W RF power) are used to reduce argon-based interferences and to increase the sensitivity of the analysis. A gas blank subtraction is performed on the data acquired at each position, based on the signal acquired 10 s prior to the ablation. Precise and accurate trace-element concentration data were obtained from offline calibration, using seven international natural and synthetic glass and carbonate reference materials: BHVO-2G, BIR-1G, GSD-1G, GSE-1G and MACS-3 (United States Geological Survey) as well as SRM 610 and 612 (National Institute of Standards and Technology). Ca is used as an internal standard for calibration of the speleothem measurements, following the assumption that the calcium carbonate in the speleothem contains 38 wt % Ca. Based on the reference materials and settings described, the reproducibility of the produced elemental concentration data was typically of the order of 5 % relative standard deviation (RSD). Limits of detection (LODs) are given in Table 2.

Table 2 Overview of limits of detection (LODs) of trace elements measured for this study using LA-ICP-MS.

3.2 Data processing

Frequency analysis is applied to study the variations in the different proxy signals, and allows for evaluating which of these proxies fluctuate seasonally. The benefit of frequency analysis for assessing seasonal cyclicity in a proxy in speleothems and other incremental climate archives was already recognized by Smith et al. (2009) and de Winter et al. (2017). Furthermore, the method can identify multi-annual trends or variability at the sub-seasonal level. Frequency analysis is performed using fast Fourier transformations (FFTs) of the isotopic and trace-element data in the distance domain. The data were detrended and padded with zeros. The power spectra were plotted as simple periodograms with frequencies in the distance domain (mm⁻¹) to allow for interpretation of seasonality in the data. The significance level (95 %) was evaluated using Monte Carlo noise simulations. The routine used operates in MATLAB^® and was based on the scripts provided in Muller and MacDonald (2000); (see Bice et al., 2012 for more detail).

Figure 4 Stable-isotope ratios and trace-element variations plotted against distance for P20, P17 and P16. Blue bars mark the DCL laminae. The left side represents the youngest layers. All stable-isotope ratios are expressed as ‰ VPDB, and trace-element concentrations are reported in parts per million (ppm). Red bars indicate years used for annual stack (Figs. 7 and 9).

An effective method to compare sub-annual variations of different proxies with each other is by resampling multiple annual cycles at a regular interval and stacking the individual cycles (Treble et al., 2003; Johnson et al., 2006; Borsato et al., 2007; de Winter et al., 2018). The advantage of this method is that the phase relations of the different proxies are preserved (Treble et al., 2003). Annual stacks are created based on moving averages to diminish the influence of low-frequency noise on the annual stacks. The number of points used for moving averages is determined as a function of the sampling resolution (i.e. three-point moving average for stable-isotope records and five-point moving average for trace-element records; see Fig. 4). Proxy records with well-constrained seasonal variation are used to define seasonal cycles. In this study, individual years were selected based on δ¹³C (minima) for stable-isotope records and Zn (maxima) for the trace-element records. Stable-isotope ratios and trace-element stacks were created separately (Fig. 3). For P16 and P17, all annual cycles were included in the stack, except for the first and the last one, since there was no guarantee that these are entirely represented in the record. For P20, only 10 years were selected from the full record to avoid the effect of multi-decadal variability (see Fig. 2) and to maintain an approach similar to that of P16 and P17. The years are indicated by the red line in Fig. 4.

Figure 5 Ranges of the stable-isotope (left) and trace-element data (right). For the stable-isotope ratios, the data mark the mean (diamonds) and the standard deviation (1σ) of the distribution. For the trace-element concentrations, the boxes represent the minimum and maximum values and the white diamonds mark the median. Numbers on top of the bars represent the percentage of the data that is below the calculated detection limit.

4 Results

The concentration range of each proxy measured in the three different intervals is shown in Fig. 5. For δ¹³C and δ¹⁸O, the mean values and ranges (minima to maxima) in P20 are significantly higher than those in P17 and P16. To illustrate the spread in the trace-element records, the median was used instead of the mean as the median is less sensitive to large concentration ranges and outliers. Al, Si, K, Mn, Rb and Th are not included in this study since >25 % of the data falls below the LOD. An exception was made in the case of Y of which only few data points are retained for P17 (81 % of the data are < LOD) and 18 % and 36 % of the data are < LOD in P20 and P16, respectively. However, Y data are discussed because of the clear seasonal signal shown in P20 and P16 (Figs. S3 and S5).

Figure 6 Periodograms (FFT) of δ¹³C, Mg, Zn and P measured in P16 to illustrate how the quality of a proxy to record the seasonal cycle can be studied. The red line represents the 95 % confidence level. δ¹³C is taken as a reference. The periodograms include two examples of proxies with a distinct peak in the seasonal frequency band of 0.8 mm⁻¹(Mg and Zn) and one proxy with no peak (P) in the seasonal frequency band. Periodograms for all periods are provided in the Supplement.

Records of stable-isotope ratios (δ¹³C and δ¹⁸O) and trace-element (Mg, P, Zn, Sr, Y, Ba, Pb and U) concentrations are plotted against distance in Fig. 4. The occurrence of darker laminae (DCL) in the samples is indicated by blue bands, clearly showing that these annual laminae are thicker in P16 and P17 (mean 1.135  m and 1.096 mm, respectively) compared to P20 (mean 0.382 mm). For all intervals, the seasonal cycles are well constrained by δ¹³C, with lower δ¹³C values occurring in DCL. Van Rampelbergh et al. (2014) present in their Fig. 4 the correspondence between stable-isotope compositions of the calcite samples taken by a micromill on a regular spatial sampling interval and a scan of the stalagmite. This correspondence is based on the information available in 2014 and 2015. However, more recent tools combined with the study of a larger portion of the stalagmite demonstrates that contrary to the affirmation of Van Rampelbergh et al., (2014, 2015), lower values in δ¹³C are found in the DCL. The mean δ¹³C is higher for P20 (−8.36 ‰) compared to P17 and P16 (−9.82 ‰ and −10.04 ‰, respectively). In addition, the amplitude of the individual cycles is larger in P20. Seasonal cycles in δ¹⁸O are much less pronounced. The most distinctive cycles are observed in P20 and some can be identified in parts of P17 and P16 (e.g. between 4 and 7 mm in P16 or between 3 and 7 mm in P17), while for other parts (e.g. between 7 and 11 mm in P16) they appear to be absent.

Figure 7 Annual stacks of the trace-element proxies. Y axis: concentrations (µg g⁻¹); x axis: sub-annual increment (1 represents January, 12 represents December). Blue shaded areas indicate the winter season. For the years used, see Fig. 4.

Seasonal variations are observed for Mg, Sr and Ba in all three intervals investigated (Fig. 4). In P17 and P16, the median concentrations of these elements are similar; 447 and 444 µg g⁻¹ for Mg, 51 and 45 µg g⁻¹ for Sr, and 36 and 33 µg g⁻¹ for Ba (Fig. 5). However, in P20 concentrations of Mg and Ba are slightly higher compared to the older intervals, i.e. 706 µg g⁻¹ for Mg and 46 µg g⁻¹ for Ba. This is also the case for Pb and U with concentrations significantly lower in P17 (0.14 and 0.05 µg g⁻¹, respectively) and P16 (0.14 and 0.07 µg g⁻¹, respectively) and a seasonal cycle that is less pronounced than in P20 (0.37 and 0.18 µg g⁻¹). In contrast, P16 has the highest median concentrations of Zn (54 µg g⁻¹) and Y (0.04 µg g⁻¹) and both elements display a well-defined seasonal covariation. Note that seasonal variability in trace elements in P16 is most pronounced in the first 9 mm, after which the amplitude of variability decreases. The same is not observed in the stable-isotope records. Although the concentration of Zn is lower in P20 and P17 (14 and 25 µg g⁻¹, respectively), the seasonal cycle is still present. Similar observations can be made for Y in P20 (0.02 µg g⁻¹). Within P16 and P17, maxima of Zn, Y, Sr and Ba mostly occur within the DCL.

Figure 6 shows an example of the FFT periodograms of δ¹³C, Mg, Zn and P in P16. Additional periodograms for the other elements in P16, P17 and P20 are included in the Supplement (Figs. S3–S5). The frequency analysis confirms the clear seasonal cyclicity of δ¹³C previously observed by Van Rampelbergh et al. (2014) (Fig. 6). The dominant frequency of δ¹³C in P16 is 0.8 mm⁻¹ (Fig. 6). This corresponds to a period of 1.25 mm, which is in good agreement with an observed mean layer thickness of 1.13 mm (Fig. S6). Growth rates based on neighbouring U-Th dates are slightly lower (0.94 mm yr⁻¹ on average; see Table 1), highlighting the benefit of combining multiple independent dating methods (e.g. layer counting and radiometric dating) in age models to increase the accuracy of dates and growth rate reconstructions. Because of its distinct seasonal cyclicity, the δ¹³C cycle is used as a reference to deduce whether or not other proxies record the seasonal cycle. Mg and Zn appear to track this seasonal cycle well as their periodograms contain peaks at 0.8 and 0.75 mm⁻¹, respectively, corresponding closely to the frequency of δ¹³C. For Zn, a broader double peak is observed with a main period of 1.18 mm and a smaller period of 1.02 mm. This double peak in the periodogram is caused by small variations in the thickness of the annual cycles around a mean thickness of 1.14 m with a lightly skewed distribution towards thinner layers (Fig. S6). The P record does not display any significant seasonal cycle (95 % confidence) (Fig. 4 and Fig. 6). For P20, visible layers are thinner (mean 0.382 mm) and also the variation in thickness is larger (RSD 28.9 %) compared to P16 and P17 (Fig. S6). This results in broader and less well defined seasonal peaks in the periodograms.

Figure 8 Cross plot showing the relationship between monthly average trace-element (Mg, Sr, Ba) concentrations and δ¹³C values in the three growth periods. Plots are labelled with the Pearson's r values indicating the strength and direction (phase relationship) of the correlation between the variables. Note the change in-phase relationship between δ¹³C and Sr and Ba, while δ¹³C and Mg remain positively correlated.

5 Discussion

5.1 Seasonal cyclicity in stable-isotope and trace-element records

A schematic overview of the observed changes in all proxies and their interpretation for the three intervals is provided in Table 1. Assessing the exact phasing of the seasonal cycles of different trace elements to δ¹³C and the visible layering remains challenging since (1) a multitude of factors control trace-element variations within speleothems and (2) stable-isotope ratios and the trace-element concentrations are not measured on the same exact sample localities. Whilst stable-isotope ratios were measured on micro-drilled powders, trace-element concentrations were measured using laser ablation. These records were later carefully aligned based on microscopic observation of sample positions, but due to differences in sample size and sampling density the exact location of the analyses may not fully coincide. An example of such a phase problem is the occurrence of an additional year in P16 in the trace-element curve compared to δ¹³C (Fig. 4, between 1 and 6 mm). Another example is the decrease in amplitude of trace-element seasonality which takes place in P16 around 9 mm from the start of the record. This change does not occur in the stable-isotope record and therefore complicates the comparison between records in the second part of the record. The fact that all trace-element records are effected and not the stable-isotope records suggests that this transition highlights a methodological issue, such as a difference in laser beam focus (e.g. due to a difference in polishing quality) during LA-ICP-MS measurements, highlighting the importance of sample preparation for high-resolution chemical analyses (e.g. LA-ICP-MS and micro-X-ray fluorescence (µXRF)). Nevertheless, δ¹³C minima generally occur in the DCL, suggesting a similar timing (and maybe control) on the visible laminae and δ¹³C cycles. Trace-element proxies show cyclicity with a similar frequency as the δ¹³C (Fig. 6). This in contrast to δ¹⁸O, which shows less clear seasonal cycles in P16 and P17.

Figure 9 Annual stacks of δ¹³C (black) and δ¹⁸O (red). Dashed lines mark the 2σ uncertainty. The x-axis intervals represent 1 year. Blue shaded areas indicate the winter season. For the years used; see Fig. 4.

5.1.1 Zinc, yttrium and lead proxies

In earlier monitoring studies carried out in the Père-Noël Cave (also part of Han-sur-Less cave system, Fig. 1), the presence of a late autumn increase in discharge was identified (Genty and Deflandre, 1998; Verheyden et al., 2008). In situ conductivity measurements indicated an increase in mineral content (ion concentrations) and/or organic matter during this autumnal increase in drip water discharge (Genty and Deflandre, 1998). Measurements of the drip water discharge above the Proserpine stalagmite show that in October and November, a doubling of the discharge volume occurs. This increased discharge is maintained until May, when a gradual decrease is observed (Van Rampelbergh et al., 2014). In the same period, Mg ∕ Ca and Sr ∕ Ca ratios in drip waters increase as a delayed response to the intense July/August rainfalls (Verheyden et al., 2008). The timing of the elevated discharge agrees with the theoretical water excess occurring above the cave (Genty and Quinif, 1996). At the same time, Verheyden et al. (2008) postulate that autumn flushing brings humic and fulvic acids that accumulate due to intense biological activity above the cave during spring and summer. The observed seasonal cycles in Zn, Y and Pb in the intervals studied are likely caused by this annual autumn–winter flushing. Variations in these trace-metal concentrations within annual speleothem layers have previously been attributed to the annual hydrological cycle. For instance, Borsato et al. (2007) linked the peak concentrations of F, P, Cu, Zn, Br, Y and Pb to the annual increase of soil infiltration during autumnal rainfall. Furthermore, it was suggested that the transport of such elements mainly occurs via natural organic matter (NOM) or other forms of colloidal material. Enrichments of these soil-derived elements within speleothems are believed to be associated with high drip water flow events (Fairchild and Treble, 2009). Studies have shown that trace metals, such as Cu, Ni, Zn, Pb, Y and rare-earth elements (REEs), are predominantly transported via complexation by NOM, of which the fraction size in the karstic waters ranges from nominally dissolved to colloidal or particulate (Hartland et al., 2012; Wynn et al., 2014). In the case of Zn and Pb, Fairchild et al. (2010) have shown that in Obir Cave (Austria) the visible and ultraviolet lamination forms during autumn and is enriched in Zn, Pb and P. According to Wynn et al. (2014), the correspondence of distinct Zn and Pb peaks with the autumnal laminae is compelling evidence for a high-flux transport of these trace metals with NOM.

However, in this study no distinct annual cycle within the P record is observed (Figs. 4 and 6). The difference between records of P and other soil-derived elements (e.g. Zn and Y; which exhibit clear seasonality) is difficult to explain. Phosphorus is considered soil derived as it originates from vegetation dieback (Treble et al., 2003; Baldini et al., 2002). Therefore, P has shown similar variations as observed in Zn, Y and Pb in previous studies (Borsato et al., 2007; Fairchild et al., 2010). In the Proserpine speleothem, no relation between P and other soil-derived trace elements is detected. This finding is in agreement with minimal seasonal variability observed in δ¹⁸O and δD values of cave water monitored in 2012–2014 (Van Rampelbergh et al., 2014), suggesting that seasonal changes in the epikarst, linked to water availability, were dominant over seasonal processes related to surface (soil) processes. If P was predominantly mobilized by humic and fulvic acids, one would expect a seasonal pattern that follows the autumn increase in discharge, in which these acids are supposedly enriched (as postulated by Verheyden et al., 2008). However, such a pattern is absent from our data.

An explanation proposed by Frisia et al. (2012) is that P is not derived from soil leaching but from other sources such as phosphate minerals present in the epikarst or microbiological activity. However, no data are available on the P concentrations in the host rock and no microbially induced petrographic features were observed in the Proserpine speleothem, making it hard to test this hypothesis.

Alternatively, the lack of seasonality in P might be explained by the occurrence of snowmelt in the earlier growth periods (P16 and P17). The delay of peak infiltration caused by snowmelt dilutes the autumn flushing effect and explains the lack of an autumn peak in the P record. Snow melt may have the so-called “piston effect” pushing older water out of the epikarst and therefore increasing the flow of water to the speleothem site, causing increased seasonality in Mg and Sr. However, there is no evidence of this delayed infiltration caused by snow melt in modern cave monitoring (Van Rampelbergh et al., 2014), while seasonality is also absent in the P record of P20. Therefore, snowmelt cannot fully explain the trace-element patterns observed in the data.

Another alternative explanation might be that the limit of detection of P in our LA-ICP-MS data is higher relative to the measured values than that of other elements of interest (e.g. Mg, Zn and Sr), causing higher analytical noise on the P record compared to the other trace-element records. The reason might be that P measurements are sensitive to interferences in a Ca-rich matrix (such as calcium carbonate in speleothems).

Because of the distinct signature of the seasonal cycle in Zn, the Zn peaks are used as tie points to create the annual stacks of other trace-element records (Fig. 7), with lower concentrations occurring during periods of lower discharge and vice versa. The much higher Zn and Y peaks in P16 compared to P17 and P20 suggest an increased seasonality effect in discharge; therefore, the accompanied annual flushing of the soil above the cave appears more intense in the early 17th century. Concentrations of Pb are significantly higher in P20 compared to the other periods (median of 0.37 µg g⁻¹ versus 0.14 µg g⁻¹ and 0.14 µg g⁻¹ in P16 and P17, respectively). An increase in Zn and Y in P20 similar to that in Pb is not observed, suggesting that the Pb-enrichment occurs at the soil level from another source. A study of Allan et al. (2015) on Pb isotope ratios in the same Proserpine stalagmite shows that the Pb concentrations are soil derived and originate from various sources of anthropogenic atmospheric pollution (coal, industrial activities, steel production and road dust). This explains well the observed higher Pb concentration in P20. Allan et al. (2015) identified increases in Pb concentration during 1945–1965 CE and 1975–1990 CE, which are in agreement with the observed higher Pb concentrations in this study between 20–18 and 13–5 mm. They also concluded that this 20th century anthropogenic pollution only affects Pb and none of the other elements used as palaeoseasonality proxy in this study.

5.1.2 Magnesium, strontium and barium proxies

Figure 7 shows that the annual stacks of Sr and Ba correlate strongly within all three intervals, evidenced by Pearson correlation coefficients (r) of 0.71, 0.97 and 0.82 for P20, P17 and P16, respectively, with p values much smaller than 0.01 (99 % confidence level; see also Fig. 8). Magnesium displays an antiphase relationship with Sr and Ba in P16 ($r=-\mathrm{0.85}$, $p=\mathrm{1.8}\times {\mathrm{10}}^{-\mathrm{7}}$), whereas in P20 this relationship is in phase (r=0.64, $p=\mathrm{7.2}\times {\mathrm{10}}^{-\mathrm{4}}$). For P17, there is no significant relationship between Mg with Sr and Ba ($r=-\mathrm{0.13}$, p=0.53).

A strong covariation of Mg with Sr and Ba, as observed in P20, has previously been attributed to reflect the presence of prior calcite precipitation (PCP) in the epikarst above the cave, caused by the occurrence of drier periods (Fairchild et al., 2000), even on a seasonal scale (Johnson et al., 2006). The presence of PCP during late summer (with high evapotranspiration above the cave), when strongly reduced drip water discharge exists above the Proserpine stalagmite, has also been evoked to explain the enriched δ¹³C of freshly deposited calcite during the cave's summer mode (Van Rampelbergh et al., 2014). Despite the difficulties of accurately correlating trace elements and stable-isotope proxies, there appears to be a distinct covariation between the P20 Mg and δ¹³C record, with maxima in Mg corresponding with maxima in δ¹³C, confirming the hypothesis of PCP control on these proxies. Like PCP, a change in the morphology of the speleothem itself, resulting in prolonged degassing, might also partly explain the correlation between δ¹³C on the one hand and Mg and Sr on the other. Both these processes are likely driven by a reduction of infiltration between the 17th and 20th century.

The antiphase correlation between Mg on the one hand and Sr and Ba on the other hand in P16 requires the involvement of other processes that dominate over PCP. A positive relationship between the Mg partition coefficient and temperature would be expected from thermodynamic considerations, and this has indeed been observed in experimental carbonate precipitation studies (Gascoyne, 1983; Rimstidt et al., 1998; Huang and Fairchild, 2001; Day and Henderson, 2013). In similar experiments, strontium partitioning into inorganic carbonate is known to remain constant with increasing temperatures but can be influenced by calcite precipitation rate (Day and Henderson, 2013). Faster precipitation of calcite causes an increased amount of lattice defects, resulting in an increased value for the partition coefficient of Sr (Pingitore, 1978; Pingitore and Eastman, 1986; Pingitore et al., 1992) and thus more Sr uptake in the calcite. Higher temperatures, combined with a decrease in drip water discharge, leading to decreased growth rates, could therefore theoretically explain the antiphase relationship of Mg and Sr (Huang et al., 2001). However, growth rates in P16 are rather high and, additionally, it has been suggested that the variations of Sr and Mg in drip water chemistry are often significantly higher than those caused by the processes mentioned above (Van Beynen et al., 1997). Roberts et al. (1998) concluded that the temperature dependence of the Mg partition coefficient could theoretically explain seasonal Mg variations but not the multi-annual trends, for which hydrological changes are likely more important. Such observations have caused the interpretation of the Mg proxy to shift from a temperature relationship to an interpretation in terms of hydrological changes such as amount of water recharge in the epikarst (Fairchild and Treble, 2009).

In this case, an alternative explanation for the P16 antiphase relation in Mg, Sr and Ba is the incongruent dissolution of dolomite (IDD of CaMg(CO₃)₂), taking place during annual periods that are characterized by enhanced water–rock interaction. The presence of dolomite within Givetian limestone deposits overlying the cave has been recognized by Verheyden et al. (2000), Fairchild et al. (2001) and Pas et al. (2016). During periods of decreased recharge, i.e. summer for the Han-sur-Lesse cave, prolonged interaction between water and rock leads to saturation of the karstic water with respect to CaCO₃. When saturation is reached, incongruent dissolution of dolomite (IDD) will start and Ca²⁺ concentration remains constant due to the precipitation of calcite (Lohmann, 1988). IDD increases the Mg ∕ Ca of the drip water (Fairchild et al., 2000) but lowers the Sr ∕ Ca and the Ba ∕ Ca ratios, because dolomite tends to have lower Sr and Ba contents with respect to calcite (Roberts et al., 1998). During winter recharge, saturation of the water in the epikarst with respect to calcite is not attained and dolomite does not dissolve. The IDD process may overwhelm the PCP signal in P16 and might be responsible for the observed antiphase relation. Since the host rock of the Han-sur-Lesse cave contains both dolomite and shale deposits (rich in trace elements such as Mg, Sr, Ba and Y), increased host rock interaction of infiltration waters may influence the trace-element concentrations of drip waters (see Fig. 1). The dominance of the interaction with dolomite or shales on the trace-element budget depends on the composition of the host rock and the local pathway of the infiltration water. This makes it very hard to separate these different processes that influence the trace-element composition of the Proserpine speleothem. The fact that host rock interaction is dominant in P16 but not in P17 and P20 might indicate that in the 16th century these periods of decreased recharge of the epikarst were more common, suggesting that a more seasonal precipitation or evapotranspiration regime was in place in the 16th century and that seasonality in water availability (and therefore epikarst recharge) was reduced in later times. Higher recharge in the 16th and 17th centuries is also evident from the thicker annual laminae in P16 and P17 compared to P20, showing that the Proserpine speleothem grew faster during these times, which points towards faster drip rates.

A third process that could cause changes in the phase relationship between Mg, Sr and Ba is a change in land use, which changes soil thickness and influences the rate of infiltration of meteoric water into the epikarst. Both Mg and Sr may have a twofold provenance: from soil and from host rock (Rutlidge et al., 2014; Belli et al., 2017). A change in land use resulting in a change in the relative contribution of soil-derived dissolved ions could therefore cause the shifts in phase relationship between Mg, Sr and Ba observed in Figs. 7 and 8. Increased seasonal cyclicity in Y and Zn, elements associated with soil activity, in P16 compared to P17 and P20 actually supports this hypothesis. The soil above the Han-sur-Lesse cave is very thin (∼25 cm), but its thickness might have varied through time. A change in vegetation cover between the 16th century and later centuries due to, for example, changes in regional climate or an increase in anthropogenic activity (e.g. forest disturbance) could have brought about such a change in soil thickness above the cave and explain the changes in trace-element patterns observed in the Proserpine speleothem.

The comparison of the annual stacks for Mg, Sr and Ba of the different intervals corroborates the idea that PCP is the main process controlling the seasonal variations of these trace elements in P17 and P20 based on the in-phase relation of Mg, Sr and Ba. Within P16, enhanced seasonality in recharge causes IDD or increased concentrations of soil-derived trace elements to dominate over PCP. This is in agreement with the lower δ¹³C values for P16, indicating a higher contribution of soil CO₂ and explaining the antiphase relation of Mg against Ba and Sr as a consequence of IDD. Somewhere between the P16 and P17 periods, a turnover in the hydrological regime of the epikarst allowed PCP to become dominant over IDD or soil contribution in the seasonal variations in the proxies. Within P17, the relationship between Mg, Sr and Ba is less clear. This could point towards a change in hydrological regime between the periods of deposition of P16 and P20, which was still underway during deposition of the calcite in section P17 (mid-17th century).

5.1.3 Uranium

In speleothems, U is thought to be mainly derived from bedrock dissolution (Bourdin et al., 2011; Jamieson et al., 2016) and to be subsequently transported by the ground water towards the speleothem (Fairchild and Baker, 2012). The partition coefficient of U is <1 for calcite (Johnson et al., 2006; Jamieson et al., 2016). This causes U to be preferentially excluded from the calcite and enriched in the remaining drip water during the process of PCP. However, in P20, where PCP is evoked as the dominant process controlling Mg, Sr and Ba seasonal variations, an antiphase relationship of U with Mg, Sr and Ba is observed (Fig. 7). Johnson et al. (2006) concluded that scavenging of U as uranyl ion (UO${}_{\mathrm{2}}^{\mathrm{2}+}$) from the drip water onto the calcite crystal surfaces during PCP has a more dominant control on seasonal U variability than the partition coefficient. Furthermore, the relatively fast growth rate of the Proserpine speleothem limits uranyl leaching from the exposed crystal surfaces in the newly formed speleothem, highlighting the dominance of epikarst processes over those taking place at the speleothem surface (Drysdale et al., 2019). Such a mechanism may explain why U is antiphase with the Mg, Sr and Ba variations, especially since the stalagmite is under continuous dripping of water.

5.2 Seasonal variations in δ¹³C and δ¹⁸O

To compare and understand the seasonal variations in δ¹³C and δ¹⁸O, annual stacks were created (Fig. 9) by resampling based on the occurrence of peaks in δ¹³C values as this proxy reflects the seasonal cycle best (Fig. 6). The minima in δ¹³C always occur in DCL in P16 and P17. In P20, this relationship is less clear, while most of the δ¹³C minima do also occur within the DCL (Fig. 4). Van Rampelbergh et al. (2014) observed seasonal changes in δ¹³C of drip water and recent calcite, with higher values occurring in summertime when drip rates are reduced. This led them to conclude that δ¹³C is mainly driven by changes in PCP. Higher δ¹³C values occur when more PCP is observed, i.e. during periods of lower recharge. The in-phase variations in Mg, Sr and Ba in P20 described above support the hypothesis of a seasonally changing degree of PCP (see Fig. 8). Seasonal variations in the amount of PCP and its effect on δ¹³C has previously been recognized in monsoon regions (Johnson et al., 2006; Ridley et al., 2015). In P16, seasonal changes in host rock interaction or soil contribution dominate the trace-element variations in Mg versus Sr and Ba over PCP. However, since the main source of carbon in Han-sur-Lesse cave waters is the vegetation cover above the cave (Genty et al., 2001), IDD is not expected to change the δ¹³C signal. For example, a case study carried out by Oster et al. (2014) showed that an increase in IDD did not affect the δ¹³C of the speleothem significantly, despite a difference of ∼0.5 ‰ in δ¹³C between the limestone and dolomite component in the host rock. Since δ¹³C is not affected by IDD, the influence of PCP on the δ¹³C remains observable. Increased degassing related to PCP increases δ¹³C in the summer season in all growth intervals. Indeed, similar as in P20, for both P17 and P16 δ¹³C minima occur within DCL, suggesting that these DCL layers in these intervals were deposited during seasonal periods of increased drip water discharge.

Observations from cave monitoring have shown that seasonal changes in cave temperature (11–15 ^∘C) are the main driver of δ¹⁸O variations in freshly deposited calcite (−7.0 ‰ to −6.2 ‰; Van Rampelbergh et al., 2014). The δ¹⁸O periodograms show that the seasonal δ¹⁸O cycle is less developed compared to δ¹³C (Fig. 6 and Figs. S3–S5). This is also expressed in the annual stacks (Fig. 9). For P20, there is a tendency towards a positive correlation of δ¹³C and δ¹⁸O but in P17 and P16 this is unclear. Although analyses of recent calcite have clearly shown that δ¹⁸O values are at least partly controlled by the cave temperature, interpretation of the seasonal δ¹⁸O changes is difficult due to the reduced seasonal cyclicity in the δ¹⁸O records compared to other proxies. However, mean δ¹⁸O values of speleothem calcite are obviously more depleted for P17 and P16 compared to P20 (Figs. 4 and 5). Temperature changes over this period could in theory explain some variation in δ¹⁸O but only if these cooler temperatures resulted in cooler cave temperatures as well and only if fractionation was in equilibrium with ambient air, which is not certain for the 17th century part of the Proserpine. In addition, this local temperature trend of ∼1 ^∘C is not visible in the long-term δ¹⁸O record of Fig. 3, suggesting that other factors have a stronger influence on the oxygen isotope fractionation in this speleothem (e.g. the δ¹⁸O of precipitation). Therefore, the hypothesis put forward here is that the lower mean δ¹⁸O values of P16 point towards an increase in winter precipitation above the cave, since Van Rampelbergh et al. (2014) has shown that winter precipitation, such as the presence of snow, above the Han-sur-Lesse cave causes a severe decrease in δ¹⁸O of the precipitation. Subsequently, this decrease is then transferred to the drip water and into the speleothem calcite.

5.3 Variability in the seasonal cycle

The observed changes of the seasonal variations in δ¹³C, Mg, Sr and Ba between P20, P17 and P16 can only be explained by a change in the processes controlling the seasonal variability in δ¹³C, Mg, Sr and Ba (Fig. 9). In the recent period, between 1960 and 2010 CE (P20), PCP is identified as the main driver for seasonal changes in Mg, Sr, Ba trace-element concentrations. This hypothesis is supported by the δ¹³C variations. In the 16th century, Mg, Sr and Ba variations suggest that IDD or soil activity, rather than PCP, dominate the seasonal signal. Fairchild and Baker (2012) defined the term transfer function to describe the quantitative relation between speleothem chemistry and changing cave environments or climate. This transfer function depends on several chemical and environmental variables, whose importance may change over time. In this case, there is a (qualitative) change in the transfer function which causes PCP to outweigh IDD and soil activity from the 16th century to modern times. This change in transfer function is likely controlled by a change in the environment around the cave (e.g. change in land use or precipitation regime) since there are no indications for drastic changes in cave morphology over the last 500 years, as interpreted from the long-term stable-isotope ratio record (Fig. 3; Van Rampelbergh et al., 2015). As mentioned above, a likely candidate for this environmental change is a change in land use (either natural or anthropogenic) between the 16th century and modern times, which could strongly influence soil activity, making soil activity a much less important factor in the trace-element budget in the 17th and 20th century compared to the 16th century.

It is known that the strength of the acting transfer function can be used as a palaeoclimate proxy. For example, Jamieson et al. (2016) demonstrated that the seasonal (anti-)correlation between δ¹³C and U ∕ Ca varies through time within a Common Era stalagmite from Belize. During drier years, reduced seasonal variability in prior aragonite precipitation causes U ∕ Ca and δ¹³C to correlate more positively compared to wetter years. This illustrates how a transfer function can be regarded as a valuable palaeoclimate proxy. In any case, a certain external forcing is necessary for a change in transfer function to take place. A prerequisite for PCP to occur is the presence of sufficient karstic voids filled with a gas phase characterized by a lower pCO₂ than that with which the infiltrating waters previously equilibrated or similarly lower pCO₂ concentrations in the cave that cause the formation of a large stalagmite above the speleothem (Fairchild and Treble, 2009). The presence of such karstic voids is dependent on the multi-annual to decadal recharge amount of the karstic aquifer. Indeed, the mean values of trace-element concentrations imply increased water availability during P16 and P17 compared to P20. More specifically, peaks in soil-derived trace-element concentrations (Zn and Y) are higher for P16, pointing towards enhanced flushing and an increased seasonality in water availability as well as an increase in soil activity. This hypothesis is supported by the observation of marked increases in drip rates in the cave during winter (Genty and Deflandre, 1998). An anthropogenic influence explains the higher concentrations of Pb in P20, but it may also partly explain the change that occurs between the 16th century and later times. Sadly, no detailed information is available about the local changes in land use in between the 16th and 17th centuries. However, it is inferred that forest cover above the Han-sur-Lesse cave system was indeed reduced between the Little Ice Age and modern times due to increased anthropogenic activity (e.g. Van Rampelbergh et al., 2015). These changes mostly occurred on the more humid slopes rather than the drier tops of the hills. Since the Salle-du-Dôme is situated under the top of the hill (formed by the anticline structure mentioned in Sect. 2.1), changes in local forest cover are relatively small, but cannot be neglected. Lastly, layers in P16 and P17 are up to 3 times thicker compared to P20 (Table 2, Figs. 3 and S6), which reflects higher growth rates. The positive relationship between water supply and growth rate has been demonstrated in the past (Baker et al., 1998; Genty and Quinif, 1996). In large and irregular shaped stalagmites, such as the Proserpine, within-layer thickness can often be quite large (Baker et al., 2008). The long-term layer thickness evolution shows a clear difference between the 17th century and present day. The significantly thinner layers during recent times clearly indicate that less water is available compared to the 17th century.

A straightforward explanation for the observed wetter cave conditions during, in particular, P16 is an increase in seasonal water excess. Gentry and Deflandre (1998) already demonstrated the strong correlation between water excess and drip rates in the Han-sur-Lesse cave system. Elevated water excess can be caused by an increase in precipitation or a decrease in temperature. Lower temperatures, especially during summer, result in a decreased evaporation of surface water. Calculations of present-day potential evapotranspiration (PET) with the Thornthwaite equation (Thornthwaite and Mather, 1957) for the period 1999–2012 show a negative water excess lasting from May to September (Fig. 2). Although the Thornthwaite and Mather (1957) method does not include vegetation effects, it is still a reliable tool to provide an idea of the effect of changes in the temperature and/or precipitation on the PET (Black, 2007). The effect of a temperature decrease during summer months on the water excess was simulated with an arbitrarily chosen 1 ^∘C temperature drop compared to the 1999–2012 mean monthly temperature. Such a temperature drop appears to have only a minor influence (Fig. 2). A hypothetical increase in total annual rainfall of 200 mm, equally spread across 12 months, has a much larger effect on the water excess (Fig. 2). However, this would decrease the length of the annual interval during which no recharge occurs (i.e. only during June–July instead of May–September) providing less suitable conditions for dolomite dissolution to occur. It must be noted that rapid channelling of excess precipitation from intense rainfall events (especially during summer) into karst voids may decrease the relevance of seasonal changes in evapotranspiration (e.g. White and White, 2013). This further stresses the importance of changes in the precipitation regime on seasonality in water excess. The land-use change proposed above would have influenced the infiltration regime. A change in forest cover since the 16th century, affecting soil thickness, would have the effect of increasing seasonality in soil-derived trace elements, especially in combination with an increase in precipitation seasonality. Therefore, two explanations could explain the seasonal-scale variability in trace-element and stable-isotope compositions observed in the Proserpine speleothem. Firstly, a stronger seasonal distribution in the amount precipitation (with more winter precipitation) could have driven an increase in host rock interaction in the 16th century, whereas today very little seasonality in the amount of rainfall is observed. Secondly, a change in soil activity due to increased forestation in the 16th century compared to modern times could have produced the strong seasonality in soil derived trace elements in Proserpine due to seasonally enhanced leaching. The presented data do not allow for a conclusive distinction between these two processes, and a combination of both (one enhancing the other) is also possible.

5.4 Implications for 17th century palaeoclimate

The majority of Common Era continental palaeoclimate reconstructions are based on tree-ring data (D'Arrigo et al., 2006), although other records, e.g. historical documents (e.g. Dobrovolny et al., 2010), ice cores (e.g. Zennaro et al., 2014) or speleothems (e.g. Baker et al., 2011; Cui et al., 2012), are used as well. Over the last decades, the consensus that changes in solar irradiance and volcanic activity are the main drivers of short-term (decadal to centennial) natural climate variability during the last millennium has been reached (e.g. Crowley, 2000; Bauer et al., 2003). Interpretations of the stable-isotope and trace-element proxies obtained on the Proserpine speleothem as well as the increased thickness of annual laminae show that a higher recharge state of the karstic aquifer characterizes the 17th century intervals compared to 1960–2010. Such an increase in recharge requires a decrease in evapotranspiration, which can result either from lower summer temperatures or higher total annual precipitation. Although it is difficult to discriminate between both the effect of a total annual precipitation increase on the recharge is expected to be higher compared to a decrease in summer temperature (Fig. 2). As mentioned above, a change in land cover (both from anthropogenic or natural causes) could also have played a role by locally affecting the infiltration of this excess precipitation. Globally dispersed regional temperature reconstructions indicate that multi-decadal warm or cold intervals, such as the Medieval Warm Period or the Little Ice Age (LIA), are not global events. Yet, a global cooling trend starting at 1580 CE is observed in the majority of the reconstructions (PAGES 2k Consortium, 2013). Several palaeoclimate reconstructions agreed upon the occurrence of a cold period around 1600 CE, with negative temperature anomalies persisting in Europe at decadal and multi-decadal scales (Ljungqvist et al., 2012; Luterbacher et al., 2016; Masson-Delmotte et al., 2013). Reconstructions of European summer temperature provided by Luterbacher et al. (2016) indicate that the coldest 11 and 51 year period since 755 CE in the area of Han-sur-Lesse cave occurred within the 17th century. These reconstructions showed a summer temperature decrease of 1–1.5 ^∘C around 1600–1650 CE. Although the 17th century has been recognized as the coldest of the past 12 centuries, hydrological climate conditions appear close to the long-term mean (Ljungqvist et al., 2016), with no significant wetting or drying trend. However, to account for the differences between the 1960–2010 interval and the 17th century observed in this study, an increase in the amount winter precipitation is needed, suggesting that climatic conditions were wetter during that time. Such a hypothesis is also supported by the depleted δ¹⁸O values in P16, indicating an increase in winter precipitation.

5.5 Implications for speleothem palaeoclimate studies

This study complements a longer speleothem stable-isotope compositional time series with shorter, higher-resolution stable-isotope and trace-element records which capture changes on a sub-annual scale. One of the major advantages of this approach is the ability to study phase relationships between proxies at the seasonal level. The seasonal cycle is the strongest cycle in Earth's climate, and therefore allows relationships between speleothem proxy records to be tested in a context that is more familiar than that of decadal-to-millennial oscillations, which are less well understood. Comparing the seasonal expression of speleothem proxies also allows for a more straightforward comparison of proxy records with cave-monitoring time series, which typically run on seasonal timescales (see Van Rampelbergh et al., 2014). The high-resolution analyses carried out within the context of this study are relatively labour intensive. Measurements at this resolution are thus likely not feasible along the full (centuries- to millennia-long) growth period of a speleothem, even if seasonal variability is consistently recorded for the entire growth duration. Instead, we propose that the application of high-resolution multi-proxy transects placed on strategic places along the growth axis of a speleothem (e.g. those parts where seasonal lamination is particularly well expressed) may be used to provide snapshots of seasonal proxy variability superimposed on longer term (decadal to millennial scale) variability in the record. This exercise teaches us that the expression of trace-element records should not be seen as a result of a constant transfer function, but rather as a complex interplay between chemical, speleological and environmental variables whose influence on the chemistry and (micro)morphology of the speleothem changes over time. Such changes in the transfer function have implications for the interpretation of the longer lower-resolution proxy records. This study therefore highlights the importance of including multiple proxies (e.g. trace-element in combination with stable-isotope, sedimentary and/or crystallographic analyses) as well as multiple sampling densities (high- and low-resolution sampling) to reliably interpret speleothem archives in terms of climate and environmental evolution.

6 Conclusions

This study of annual trace-element and stable-isotope (δ¹³C and δ¹⁸O) variations over three different time intervals of the annually laminated Proserpine stalagmite from the Han-sur-Lesse cave (Belgium) shows that seasonal changes in Mg, Sr and Ba during the recent period (1960–2010) suggest a strong effect of prior calcite precipitation, caused by lower water availability during summer. In the 17th century (1600 CE ± 30 and 1640 CE ± 30), however, Mg is in antiphase with Sr and Ba. This implies that another process overwrites the PCP-dominated seasonal cycle in these trace elements. A varying degree of incongruent dolomite dissolution, with more dissolution occurring during summer when water residence times in the epikarst are longer, or a more dominant influence of soil activity on the trace-element budget due to changes in land use are plausible hypotheses. The transfer function governing the trace-element concentrations in the Proserpine speleothem is driven over the last centuries by changing contributions of PCP, dolomite dissolution and soil leaching. Stable-isotope ratios (δ¹³C and δ¹⁸O), soil-derived trace-element concentrations (Zn, Y and Pb) and speleothem morphology indicate that the multi-annual recharge of the epikarst was higher in the 17th century. The change in the response of Mg, Sr and Ba in the Proserpine speleothem to environmental changes was identified to be driven by climate or by changes in land use. It may result from an increase in recharge caused by a combination of lower summer temperatures and an increase in the amount of winter precipitation in the 17th century for the Han-sur-Lesse cave region or by a change in the vegetation cover between the 16th and 17th century, which reduced the importance of soil processes on trace-element compositions in the cave's drip waters. The effect of an increase in winter annual precipitation on the recharge is expected to be larger compared to a decrease in summer temperature. The data obtained in this study therefore point towards a stronger seasonal cycle in cave hydrology during the 17th century.

This high-resolution multi-proxy study provides a good example of how the relative importance of different processes on trace-element concentrations in speleothem calcite can change over time. This observation has implications for future speleothem-based palaeoclimate reconstructions, since transfer functions for specific cave sites, determined by cave monitoring, are often assumed to remain constant when no drastic changes in the cave environment have occurred. As the change in trace-element proxy transfer function observed in this study is driven by environmental change, this change by itself can serve as a valuable palaeoclimate proxy.