Comparison of paleobotanical and biomarker records of mountain peatland and forest ecosystem dynamics over the last 2600 years in central Germany

As peatlands are a major terrestrial sink in the global carbon cycle, gaining an understanding of their development and changes throughout time is essential in order to predict their future carbon budget and potentially mitigate the adverse outcomes of climate change. With this aim to understand peat development, many studies have investigated the paleoecological dynamics by analyzing various proxies, including pollen, macrofossil, elemental, and biomarker analyses. However, as each of these proxies is known to have its own benefits and limitations, examining them in parallel allows for a deeper understanding of these paleoecological dynamics at the peatland and a systematic comparison of the power of these individual proxies. In this study, we therefore analyzed peat cores from a peatland in Germany (Beerberg, Thuringia) to (a) characterize the vegetation dynamics over the course of the peatland development during the late Holocene and (b) evaluate to what extent the inclusion of multiple proxies, specifically pollen, plant macrofossils, and biomarkers, contributes to a deeper understanding of those dynamics and interaction among factors. We found that, despite a major shift in the regional forest composition from primarily beech to spruce as well as many indicators of human impact in the region, the local plant population in the Beerberg area remained stable over time following the initial phase of peatland development up until the last couple of centuries. Therefore, little variation could be derived from the paleobotanical data alone. The combination of pollen and macrofossil analyses with the elemental and biomarker analyses enabled further understanding of the site development as these proxies added valuable additional information, including the occurrence of climatic variations, such as the Little Ice Age, and more recent disturbances, such as drainage.

1 Introduction

Peatlands have been well established as an important sink in the terrestrial carbon cycle, containing about 25 % (600 Gt) of the global soil carbon stock (Yu et al., 2010), though only comprising 3 % of the global land area (Xu et al., 2018). The same characteristics that make peatlands an effective carbon sink, e.g., slow degradation rates and high accumulation rates of organic matter, also make them an excellent archive for paleoecological records at high temporal resolution (Barber et al., 1994; Chambers et al., 2012). Compared with other environmental archives, the use of peat sequences, particularly ombrotrophic bogs, for investigating paleoenvironmental changes is advantageous as they are generally more accessible; contain material readily dated using radiocarbon, resulting in chronologies with lower uncertainties; and are primarily influenced by atmospheric inputs only, thereby recording climatic information along with their ecologic response (Chambers et al., 2012). Furthermore, peat sequences provide some of the few well-resolved records of paleoenvironmental changes during the Holocene in central Europe, aside from lake sediments (e.g., Schwark et al., 2002; Hepp et al., 2019).

Through extensive studies across ecosystems, including peatlands, multiple proxies have been developed and applied to characterize paleovegetation and paleoclimatic conditions, including pollen (e.g., Speranza et al., 2000), macrofossils (e.g., Tuittila et al., 2007), and biomarkers (e.g., Ficken et al., 1998a; Xie et al., 2004). Although these proxies have all been reliably used in paleoenvironmental reconstructions, each has benefits and drawbacks in their application.

Pollen and plant macrofossils are commonly used paleobotanical analyses to (a) reconstruct local and regional plant succession of peatlands and their surroundings, (b) assess human or climate impact on plant communities and understand their response or resilience, and (c) define pristine plant populations as reference conditions in their restoration (Speranza et al., 2000; Gałka et al., 2022b). Pollen records extracted from peatlands have a well-established history of use for reconstructing paleovegetation, beginning with Von Post's study of arboreal pollen in Swedish bogs (1916). Despite the benefits of using pollen to reconstruct vegetation, there can be several limitations. Pollen can be transported over large distances, meaning that the local vegetation signal is overlain by that of the region (Farrimond and Flanagan, 1996). Additionally, some plants produce more pollen than others, and some pollen is better preserved, leading to over-representation of certain taxa in the results (Birks and Birks, 2000). These drawbacks can be mitigated by including macrofossil analysis alongside pollen (Birks and Birks, 2000). Due to waterlogged conditions and low pH, the remains of mosses, graminoids, and dwarf shrubs are well preserved in peatlands, allowing for the identification of macrofossils to reconstruct the development and succession of in situ, peat-forming vegetation within peatlands (Chambers et al., 2012). Macrofossils often provide better taxonomic precision than pollen (Birks and Birks, 2000), which is important when reconstructing the potential of individual plants to aid in carbon sequestration (Sim et al., 2021). However, macrofossils can also be difficult to identify in highly humidified peat (Naafs et al., 2019). In such cases, or if changes are more subtle, plant biomarkers can be a valuable addition to studies of peat cores.

Plant biomarkers derived from leaf waxes have been used extensively in paleoecological investigations with a long history of research throughout the past century (e.g., Chibnall et al., 1934). The most commonly used of these markers are the straight-chain lipids: n-alkanes, n-alkanols, and n-fatty acids (Jansen and Wiesenberg, 2017). Previous studies of peatlands have used these biomarkers as proxies for paleoclimate and vegetation (e.g., Ficken et al., 1998a; Xie et al., 2004). As peat consists almost entirely of organic matter, there are large concentrations of a range of biomarkers. While biomarkers may also be deposited as aerosols to a small extent, the vast majority derive directly from the parent peat vegetation, and even as the peat becomes more humified, biomarkers have the potential to persist (Naafs et al., 2019). Aliphatic compounds, including n-alkanes, are particularly resistant to decomposition and become residually enriched (Biester et al., 2014). In peatlands particularly, n-alkanes have been used to differentiate between vegetation types, as Sphagnum mosses typically produce shorter carbon chains (23–25), while longer carbon chains (27–33) are commonly produced by vascular plants (Chambers et al., 2012; Pancost et al., 2002; Baas et al., 2000). Indices comparing relative abundances of biomarker homologues, including the carbon preference index (CPI), can also be used to indicate degradation of organic matter as well as climate changes that are more favorable for high microbial activity and enhanced degradation (Chambers et al., 2012). Yet, the use of biomarkers and biomarker-derived indices and ratios to reconstruct vegetation and other environmental conditions is complicated by the fact that there are multiple sources of biomarkers, and species-specific compositions are rare (Jansen and Wiesenberg, 2017). It is more straightforward to differentiate between vegetation types, but even this can be complicated due to factors that influence biomarker distribution patterns, such as moisture, humidity, and the stage of plant and/or leaf development (Jansen and Wiesenberg, 2017).

Previous studies have successfully combined pollen and lipid biomarker analyses to reconstruct regional paleoenvironmental conditions (e.g., Farrimond and Flanagan, 1996; Schwark et al., 2002), including in peatlands (e.g., Zhou et al., 2005). Additionally, a growing number of studies have included both pollen and macrofossil analyses along with biomarker measurements for a more robust interpretation of peat archives (e.g., Ronkainen et al., 2015; Balascio et al., 2020). Considering elemental analyses of carbon, nitrogen, and their stable isotopes also allows for more nuanced interpretations as the $\mathrm{C}/\mathrm{N}$ ratio can be used as an indicator of decomposition (Kuhry and Vitt, 1996; Biester et al., 2014) and δ¹³C of decomposition and hydrological changes (Loisel et al., 2010; Biester et al., 2014).

In this study, we aimed to evaluate these proxies individually and in combination while characterizing the paleovegetation dynamics of the Beerberg peatland, an ombrotrophic mountain peat bog in the Thuringian Forest (Germany). The Beerberg site was chosen as there are few paleoenvironmental archives of the late Holocene in central Germany and more specifically in the Thuringian Forest (e.g., Breitenbach et al., 2019), leading to a demand for a better understanding of such scarce paleoenvironmental records (Githumbi et al., 2022). Local reconstructions of Thuringia are especially valuable as the area lies near the current boundary between the maritime Cfb (temperate, no dry season, warm summer) and continental Dfb (continental, no dry season, warm summer) climatic zones in the Köppen–Geiger classification (Peel et al., 2007). This boundary shifts spatially through time, and increased dated and high-resolution reconstructions are needed to fine-tune knowledge of the past evolution of the climate boundaries in central Europe (Breitenbach et al., 2019). Previous pollen analyses have been completed at the site by Jahn (1930) and Lange (1967) but were more limited in scope as Jahn only considered arboreal pollen, and neither study included radiocarbon dating, so there is no reliable chronology for the Beerberg peatland. Based on the previous point, the aim of the present study was to assess (a) the paleovegetation dynamics over time in the Thuringian Forest over the last ca. 2600 years, an area currently lacking such continuous records, and (b) how adding n-alkanols and n-fatty acids to a multi-proxy analysis helps to obtain a more detailed reconstruction of past environmental conditions and peatland development.

Figure 1 Map showing the sampling spot in the Beerberg peatland (area within continuous white boundary), with the contour lines and inset showing the location in continental Europe.

2 Materials and methods

2.1 Study area and sampling

The Beerberg peatland (50^∘39^′32^′ N, 10^∘44^′36^′′ E; 983 m) is a raised ombrotrophic peat bog at the summit of the Großer Beerberg in the Thuringian Forest (Germany) and is part of the Vessertal–Thuringian Forest Biosphere Reserve (Fig. 1). The underlying bedrock is rhyolite (Lützner et al., 2012). Annual precipitation is estimated to be 1300 mm (Görner et al., 1984), and the mean annual temperature is 4 ^∘C (Jeschke and Paulson, 2000). While the Beerberg peatland has been under nature protection status since 1939, previously some of the peat was used as fuel for glass manufacturing. Consequently, the peatland was partially drained and reforested with spruce trees in the 19th and 20th centuries. Later, in particular after 1990, drainage trenches were filled and spruce trees were removed to prevent further drying of the peat (Thüringer Landesanstalt für Umwelt und Geologie, 2002). Current vegetation includes Sphagnum mosses, such as S. fuscum, S. magellanicum (s.l.), S. angustifolium, and S. capillifolium; tree species Picea abies, Pinus sylvestris, Betula pendula, and Betula pubescens; and other plants common in partially disturbed bogs, such as Calluna vulgaris, Eriophorum vaginatum, and Polytrichum strictum.

Sampling was completed in October 2019. Two Russian peat corers (5 cm diameter, Eijkelkamp, Giesbeck, the Netherlands; 7 cm diameter, self-made) were used to alternately core two small hummocks within an approximately 20 cm distance to a depth of 340 cm. Overlapping core sections were taken from 90–100, 175–190, and 290–320 cm. The overlapping sections were correlated by sampling depth and stratigraphy, and for the elemental and biomarker analyses, the overlapping samples were treated as replicates, and the averages of the results were calculated.

2.2 Elemental analysis

For elemental analysis, samples were taken at 5 cm intervals, with a few exceptions due to visibly distinct layers in the peat at 10–12, 170–172, 270–272, 325–327, 327–328, and 337.5–340 cm. The samples were freeze-dried to a constant weight and subsequently milled to a fine powder using a horizontal ball mill (MM 400, Retsch). Approximately 1 mg of material per sample was used for the elemental analysis. Carbon and nitrogen concentrations (C, N), as well as stable C isotope values (δ¹³C), were measured using an elemental analyzer isotope ratio mass spectrometer (EA-IRMS; FLASH 2000 HT Plus, linked by ConFlo IV to DELTA V Plus IRMS; Thermo Fisher Scientific). Calibration was carried out using caffeine (Merck, Germany) and a soil reference material from a Chernozem (Harsum, Germany; see black carbon reference materials, https://www.geo.uzh.ch/en/units/2b/Services/BC-material/Environmental-matrices.html, last access: 20 March 2023). At least two analytical replicates were measured for all samples.

2.3 Radiocarbon dating

Handpicked plant remains (Table A1) were selected from 12 depths (7.5, 16.5, 34.5, 54.5, 69.5, 124.5, 174.5, 258.5, 278.5, 293.5, 314.5, and 334–336 cm). These were cleaned by an acid–alkali–acid treatment and combusted at 900 ^∘C to produce CO₂, which was reduced to graphite. The carbon isotope composition was measured by accelerator mass spectrometry (AMS) at the Institute of Ion Beam Physics at the Swiss Federal Institute of Technology (Zurich, Switzerland) using the 0.2 MV MICADAS facility.

2.4 Macrofossil analysis

Plant macrofossils were analyzed at 4 cm resolution (n=70, 1 cm thick peat slices) using samples with a volume of approximately 8 cm³. The samples were washed and sieved under a warm-water current over 0.20 mm mesh screens. The percentage of individual fossils of vascular plants and brown mosses was estimated, and the fossil carpological remains and vegetative fragments (leaves, rootlets, epidermis) were identified using identification keys (Smith, 2004; Mauquoy and Van Geel, 2007) and compared with recently collected specimens. S. capillifolium and Sphagnum rubellum were grouped together (as S. capillifolium/rubellum) due to the difficulty in differentiating them in fossil state. Both species are typical ombrotrophic mosses that occur together in relatively dry hummocks or lawns (Hölzer, 2010), with S. rubellum preferring moister habitats (Hölzer, 2010). Similarly, Sphagnum medium and Sphagnum divinum (in the past assigned to S. magellanicum; see Laine et al., 2018) are expressed as S. medium/divinum because it was impossible to distinguish them based on morphological features in fossil state.

2.5 Pollen analysis

For pollen analysis, the core was sampled at 2.5, 6.5, 9.5, and then 5 cm intervals starting at a depth of 12.5 cm, ending with a total of 69 samples. At each depth, a volume of 2 cm³ was extracted and subjected to standard laboratory procedures for pollen analyses (Berglund and Ralska-Jasiewiczowa, 1986). Samples were treated with 10 % hydrochloric acid (HCl) to dissolve carbonates, heated in 10 % potassium hydroxide (KOH) to remove humic compounds, and finally soaked in 40 % hydrofluoric acid (HF) for at least 24 h to remove the mineral fraction. One Lycopodium tablet (10 679 spores; produced by Lund University) was added to the samples (Stockmarr, 1971). Sample slides were analyzed using an Eclipse 50i upright microscope and counted to a sum of at least 500 arboreal pollen (AP) grains. However, due to the high peat accumulation rate, this sum was not achieved in 33 samples, including 7 samples in which a 100 AP sum was not achieved. Pollen taxa were identified using atlases (Beug, 1961; Moore et al., 1991) and the reference grains owned by the Institute of Geoecology and Geoinformation, Adam Mickiewicz University, Poznań. Selected non-pollen palynomorphs (NPPs), such as fungi and algae, and microscopic charcoal particles (size fractions: 0.01–0.1 mm; > 0.1 mm) were also counted. Microscopic charcoal particles were counted until their number summed with simultaneously counted Lycopodium spores reached 200 (Finsinger and Tinner, 2005; Tinner and Hu, 2003). Palynological indicators of human impact were organized according to Behre (1981) and Gaillard (2013).

2.6 Biomarker analysis

For biomarker analysis, the same samples were used as for the elemental analysis. Soxhlet extraction was performed using 468 mg to 758 mg of the milled peat samples, as described by Wiesenberg and Gocke (2017). Briefly, total lipids were extracted over approximately 30 h with dichloromethane(DCM) : methanol(MeOH) (93:7, $v/v$). These extracts were then separated sequentially into three fractions containing, respectively, neutral components including n-alkanes and n-alkanols, n-fatty acids, and polar and high-molecular-weight compounds. For the separation, a glass column containing silica 60 + 5 % potassium hydroxide (KOH), 63–200 µm, was used along with the solvents DCM, DCM : formic acid (99:1, $v/v$), and DCM : MeOH (1:1, $v/v$). The neutral fraction was further separated into aliphatic, aromatic, and heterocompound fractions. For this, a pasteur pipette containing activated silica gel (100 Å, 70–230 mesh, dried for at least 8 h at 110 ^∘C) was used along with the solvents n-hexane, n-hexane : DCM (1:1, $v/v$), and DCM : MeOH (93:7, $v/v$).

Prior to measurement, an internal standard of deuterated n-eicosanoic acid (D₃₉C₂₀; Cambridge Isotope Laboratories, Inc.) was added to the fatty acid fraction, which was then methylated using a boron trifluoride–methanol solution (CAS no. 373-57-9, Sigma-Aldrich). Additionally, an internal standard of deuterated n-octadecanol (D₃₇C₁₈; Cambridge Isotope Laboratories, Inc.) was added to the heterocompound fraction, containing n-alkanols, which was then silylated using N,O-bis(trimethylsilyl)acetamide (BSA) (CAS no. 10416-59-8, Sigma-Aldrich). An internal standard of deuterated n-tetracosane (D₅₀C₂₄; Cambridge Isotope Laboratories, Inc.) was added to the aliphatic fraction before analysis. The n-alkanes, n-alkanols, and n-fatty acids were quantified on a gas chromatograph (GC; Agilent 7890B) equipped with a multimode inlet and a flame ionization detector (FID). Compound identification was performed on the Agilent 6890N GC equipped with a split–splitless injector coupled to the Agilent 5973 mass selective (MS) detector. Both instruments were equipped with a DB-5MS column (50 m × 0.2 mm × 0.33 µm) and 1.5 m deactivated pre-column, with helium as the carrier gas (1 mL min⁻¹). The GC oven temperature for n-alkanes was held at 70 ^∘C for 4 min and increased to 320 ^∘C at a rate of 5 ${}^{\circ }\mathrm{C}{\min }^{-\mathrm{1}}$ held for 50 min. For n-fatty acids and n-alkanols, the temperature was held at 50 ^∘C for 4 min, then increased to 150 ^∘C at a rate of 4 ${}^{\circ }\mathrm{C}{\min }^{-\mathrm{1}}$, and finally increased to 320 ^∘C at 3 ${}^{\circ }\mathrm{C}{\min }^{-\mathrm{1}}$ held for 40 min. The samples (1 µl) were always injected in splitless mode. The gas chromatograph–mass spectrometer (GC–MS) was operated in electron ionization mode at 70 eV and scanned from $m/z$ 50–550. Individual compounds were identified through comparison of mass spectra with those of external standards and from the NIST and Wiley mass spectra library.

2.7 Data processing and analysis

2.7.1 Data entry, processing, and analysis

Data were entered, organized, and screened in Microsoft Excel. Subsequent data processing and analysis were conducted in R version 4.0.4 (R Core Team, 2023). The data were combined in five tables, one each for the radiocarbon dating, elemental analysis, plant macrofossils, pollen, and biomarker composition. In each table, each row represents one depth and each column one parameter. All data are available at Pangaea via https://doi.org/10.1594/PANGAEA.961142.

2.7.2 Elemental analysis

The mean and standard error were calculated for the analytical replicates used in the elemental analysis. To be consistent across all of the data interpretation and to determine the locations of potential climate and vegetation changes in the core, we applied a hierarchical clustering analysis most commonly performed on vegetation abundance data. We used the vegan (Oksanen et al., 2022) and rioja (Juggins, 2022) R packages. The Euclidean distance was determined using the concentrations of C and N, the $\mathrm{C}/\mathrm{N}$ ratio, and the C stable isotope values in the dist function. A constrained hierarchical clustering approach (CONISS; Grimm, 1987) was performed on the dissimilarity using the chclust function, with clusters constrained by depth. The number of zones was determined using the broken-stick model (MacArthur, 1957; Bennett, 1996) with the bstick function.

2.7.3 Radiocarbon dating

An age–depth model was developed using the Bacon model (Blaauw and Christen, 2011), as implemented in the rbacon package (Blaauw et al., 2021). The uppermost radiocarbon date was excluded from the age–depth model, as it was from a sample within 10 cm of the surface, and the yielded date of −552 yr BP ± 23 is too young and stratigraphically inconsistent with the rest of the core and represents a postbomb peak of atmospheric ¹⁴C, following nuclear weapon testing in the 1950s and 1960s (Table A1). As a constraint, an estimated surface sample age of −70 yr BP ± 5 was added. IntCal20 was used as the calibration curve (Reimer et al., 2020), and NH1 was used as the calibration curve for postbomb dates (Hua et al., 2013). All of the dates referred to in the following sections are the mean values returned by the age–depth model, which were calculated at 1 cm resolution.

2.7.4 Macrofossil

A CONISS cluster analysis (Grimm, 1987) was performed on the macrofossil abundance data. In contrast with the elemental data, the Bray–Curtis dissimilarity was calculated as this is a better measure for purely compositional data. This was completed using the vegdist, chclust, and bstick functions described in the previous section.

2.7.5 Pollen and microcharcoal

The microscopic charcoal accumulation rate (CHAR_micro; unit: $\mathrm{particles}{\mathrm{cm}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) was calculated as follows:

$$\begin{array}{}\text{(1)}& {\mathrm{CHAR}}_{\mathrm{micro}}={\mathrm{CHAC}}_{\mathrm{micro}}\cdot {\mathrm{AR}}_{\mathrm{deposits}},\end{array}$$

in which CHAC_micro is the concentration of the microscopic charcoal particles (unit: particles cm⁻³), and AR_deposits is the peat or sediment accumulation rate (unit: cm yr⁻¹) (Davis and Deevey, 1964).

Pollen percentages were calculated as taxon percentages with

$$\begin{array}{}\text{(2)}& \mathrm{\text{taxon percentages}}=(\text{number of taxon grains}/\mathrm{TPS})\times \mathrm{100}\mathit{\%},\end{array}$$

where TPS indicates the total pollen sum, including the AP and non-arboreal pollen (NAP) taxa and excluding the local taxa (i.e., aquatic, wetland, and spore-producing).

TILIA software was used to plot the diagram presenting the results of the palynological analysis (Grimm, 1993). We performed a CONISS cluster analysis as described in the previous section using the absolute pollen counts. Only taxa present in at least three samples and that reached at least 3 % relative abundance in one sample were included. Non-pollen palynomorphs and coprophilous fungi were not included in the analysis.

2.7.6 Biomarker

Biomarker amounts are reported as absolute concentrations in micrograms per gram. The total lipid extract (TLE) was calculated in milligrams per gram. The carbon preference index (CPI) (Marzi et al., 1993) and average chain length (ACL) (Poynter et al., 1989) were also calculated for the n-alkanes using the following equations:

$$\begin{array}{}\text{(3)}& {\displaystyle }\mathrm{CPI}={\displaystyle \frac{\left({\sum }_{i=n}^m{\mathrm{C}}_{\mathrm{2}i+\mathrm{1}}+{\sum }_{i=n+\mathrm{1}}^{m+\mathrm{1}}{\mathrm{C}}_{\mathrm{2}i+\mathrm{1}}\right)}{\mathrm{2}\left({\sum }_{i=n+\mathrm{1}}^{m+\mathrm{1}}{\mathrm{C}}_{\mathrm{2}i}\right)}},\text{(4)}& {\displaystyle }\mathrm{ACL}={\displaystyle \frac{{\sum }_{i=n}^m(\mathrm{2}i+\mathrm{1})\cdot {\mathrm{C}}_{\mathrm{2}i+\mathrm{1}}}{{\sum }_{i=n}^m{\mathrm{C}}_{\mathrm{2}i+\mathrm{1}}}},\end{array}$$

where C_x is the concentration of each lipid containing x carbon atoms; n and m are the chain lengths of the starting and ending lipids divided by 2, respectively (note that both 2n and 2m should be even numbers). For the n-alkanes, m is 11 and n is 15.

For the n-alkanols and n-fatty acids, the formulas were slightly adjusted to

$$\begin{array}{}\text{(5)}& {\displaystyle }{\mathrm{CPI}}_{\mathrm{ALC}/\mathrm{FA}}={\displaystyle \frac{\left({\sum }_{i=n}^m{\mathrm{C}}_{\mathrm{2}i}+{\sum }_{i=n+\mathrm{1}}^{m+\mathrm{1}}{\mathrm{C}}_{\mathrm{2}i}\right)}{\mathrm{2}\left({\sum }_{i=n+\mathrm{1}}^{m+\mathrm{1}}{\mathrm{C}}_{\mathrm{2}i+\mathrm{1}}\right)}},\text{(6)}& {\displaystyle }{\mathrm{ACL}}_{\mathrm{ALC}/\mathrm{FA}}={\displaystyle \frac{{\sum }_{i=n}^m\left(\mathrm{2}i\right)\cdot {\mathrm{C}}_{\mathrm{2}i}}{{\sum }_{i=n}^m{\mathrm{C}}_{\mathrm{2}i}}}.\end{array}$$

For the n-alkanols, m is 10 and n is 14. For the n-fatty acids, m is 10 and n is 16.

A number of indicative ratios have been developed to help interpret biomarker data, particularly for n-alkanes. We applied many of these to our data to see which was the most applicable in the Beerberg setting. The P_aq (Ficken et al., 2000) and P_wax (Zheng et al., 2007) proxies have previously been used in sediments to differentiate between aquatic macrophyte and terrestrial plant input and in peatlands to infer past water levels with high P_aq values associated with a higher water table and high P_wax values associated with a lower water table (e.g., Zheng et al., 2007; Zhou et al., 2005; Nichols et al., 2006; Andersson et al., 2011).

$$\begin{array}{}\text{(7)}& {\displaystyle }{P}_{\text{aq}}={\displaystyle \frac{{\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{25}}}{{\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}}+{\mathrm{C}}_{\mathrm{31}}}},\text{(8)}& {\displaystyle }{P}_{\text{wax}}={\displaystyle \frac{{\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{29}}+{\mathrm{C}}_{\mathrm{31}}}{{\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{29}}+{\mathrm{C}}_{\mathrm{31}}}}\end{array}$$

Andersson et al. (2011) found that P_aq and P_wax could be misleading if S. fuscum and Betula species were present in the peat due to their relatively high abundances of the C₂₃ homologue. They derived the ratio of ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ to correct for these inputs, and as both S. fuscum and Betula species were present at the Beerberg site, this ratio was included. The ratio of n-alkanes (${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$) was also calculated, as this has been used in peatland settings before to determine shifts in Sphagnum species (McClymont et al., 2008). Furthermore, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$ and ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$ were calculated, and other studies have used these to distinguish Sphagnum species in peat settings (e.g., Ronkainen et al., 2013). Far fewer ratios have been identified for n-alkanols and n-fatty acids. Nevertheless, we tested ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$, $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ (Zheng et al., 2011), and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ for the n-alkanols and ${\mathrm{C}}_{\mathrm{15}}/({\mathrm{C}}_{\mathrm{15}}+{\mathrm{C}}_{\mathrm{16}})$ (Zheng et al., 2007) and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ for the n-fatty acids.

A CONISS cluster analysis (Grimm, 1987) was performed on the absolute concentrations of the measured homologues of the n-alkanes (C₁₉–C₃₃), n-alkanols (C₁₆–C₂₈), and n-fatty acids (C₁₄–C₃₂) to determine whether different phases could be distinguished using only the biomarker data. The same procedure was followed as previously described in the macrofossil and pollen sections.

Figure 2 Age–depth model of the Beerberg core. For precise radiocarbon dates, see Table A1. The age–depth model is shown in the main graph (d). The distribution of the calibrated ¹⁴C dates is shown in blue, and the estimated surface age of −70 cal yr BP is shown in green. The wider the distribution, the less precise the dates. The dashed red curve shows the mean ages derived from the model, and the dashed grey curves represent the 95 % confidence intervals. The top three graphs in the figure show the fit of the Markov chain Monte Carlo (MCMC) iterations (a), the accumulation rate (yr cm⁻¹) with prior (green) and posterior (grey) distributions (b), and the prior (green) and posterior (grey) distributions for the memory (c) or how much the accumulation rate is able to change from one depth to the next. For more information on the Bacon model, see Blaauw and Christen (2011).

3 Results

3.1 Radiocarbon dating

The age–depth model that best matched the data in this study used a mean accumulation rate of 1 mm yr⁻¹. The resulting curve is visualized in Fig. 2. There are three clear phases with distinct accumulation rates in the model: 0.66 mm yr⁻¹ from 340–293.5 cm (2528–1826 cal yr BP), 1.99 mm yr⁻¹ from 293.5–124.5 cm (1826–978 cal yr BP), and 1.27 mm yr⁻¹ from 124.5–0 cm (978–present). Hence, the changes in accumulation rates correspond with the dated samples at 293.5 and 124.5 cm.

Figure 3 The carbon (C) and nitrogen (N) concentrations, $\mathrm{C}/\mathrm{N}$ ratio, and stable carbon isotope composition (δ¹³C) plotted against depth (cm) and modeled age (cal yr BP), including the resulting dendrogram of the CONISS analysis and phase boundaries indicated at 2113 cal yr BP (310 cm), 1569 cal yr BP (250 cm), 1151 cal yr BP (160 cm), and 809 cal yr BP (105 cm). An image of the sampled core is inset on the left. Two analytical replicates of each sample were measured. Error bars showing standard error are present for the C and N concentrations, $\mathrm{C}/\mathrm{N}$ ratio, and δ¹³C values.

3.2 Elemental analysis

The results of the elemental analysis are shown in Fig. 3. The CONISS analysis based on Euclidean distance resulted in five zones: 2528–2113 cal yr BP (Phase I-G), 2113–1569 cal yr BP (Phase II-G), 1569–1151 cal yr BP (Phase III-G), 1151–809 cal yr BP (Phase IV-G), and 809 cal yr BP–present (Phase V-G). The transitions between phases appear to be the most correlated with changes in the δ¹³C values. The carbon concentration (C) ranged from 41.1 %–51.1 %, with one exceptionally low value of 21.1 % in the basal sample of Phase I-G (340–310 cm). The nitrogen concentration (N) also had a sharp peak of 2.5 % in Phase I-G at 2368 cal yr BP (327 cm). Both the $\mathrm{C}/\mathrm{N}$ ratio and the δ¹³C values increased during Phase I-G to 64.3 and −25.8, respectively. This increasing trend continued into Phase II-G for both (310–250 cm), with the $\mathrm{C}/\mathrm{N}$ ratio averaging 90.2 and δ¹³C averaging −25.4. N values decreased, while C stayed stable. During Phase III-G (250–160 cm), the $\mathrm{C}/\mathrm{N}$ ratio fluctuated rapidly between 90.1–197.7, likely due to relatively small shifts in N. This phase contained the highest $\mathrm{C}/\mathrm{N}$ values with low N as well as less negative δ¹³C values, ranging from −23.5 to −25.2. In Phase IV-G (160–105 cm), $\mathrm{C}/\mathrm{N}$ and δ¹³C began to decrease, ranging from 95.1–159.8 and −24.1 to −25.6, respectively. In Phase V-G (105–0 cm), the $\mathrm{C}/\mathrm{N}$ ratio continued to decrease as N values increased, reaching a peak of 2 % at 102 cal yr BP (20 cm). While the δ¹³C value of the topmost sample at −29.3‰ was particularly low, there was also an interval of increasing values from 593 cal yr BP (80 cm) to 345 cal yr BP (50 cm), reaching a local maximum of −24.5.

Figure 4 Plant macrofossil diagram presenting local vegetation development in the Beerberg peatland. Taxa with percentages are estimated volume percentages; the others are absolute counts (note scale differences on the x axes). The resulting dendrogram of the CONISS analysis is on the right with phase boundaries indicated at 2251 cal yr BP (318.5 cm), 1671 cal yr BP (273.5 cm), and 64 cal yr BP (15.5 cm).

3.3 Macrofossil analysis

Based on the macrofossil analysis, the primary peat-forming species were S. fuscum and S. medium/divinum (Fig. 4). S. fuscum was dominant over most of the core, from 2086–106 cal yr BP (308.5–20.5 cm), with S. medium/divinum taking over in more recent layers (47 to −48 cal yr BP; 13.5–3.5 cm). Additionally, E. vaginatum was an important species with two major periods from 64 to −4 cal yr BP (15.5–7.5 cm) and 2528–2251 cal yr BP (340–318.5 cm).

The CONISS analysis based on the Bray–Curtis dissimilarity resulted in four significant phases (Fig. 4). Phase I-M (2528–2251 cal yr BP) is characterized by a high abundance of E. vaginatum as well as a relatively large amount of macrocharcoal (here charred wood pieces), especially at a depth of 328.5 cm (2388 cal yr BP). In Phase II-M (2251–1671 cal yr BP), S. fuscum was dominant. In Phase III-M (1671–64 cal yr BP), S. fuscum remained dominant, but there was also a small but relatively steady presence of E. vaginatum and Ericaceae rootlets. In Phase IV-M (64 cal yr BP–present), S. medium/divinum replaced S. fuscum, and there was an increase in Ericaceae rootlets as well as E. vaginatum.

Figure 5 Diagram of selected pollen taxa, NPPs, and microscopic charcoal influx (CHAR_micro), including the resulting dendrogram of the CONISS analysis and phase boundaries indicated at 1816 cal yr BP (292.5 cm), 1092 cal yr BP (147.5 cm), and 366 cal yr BP (52.5 cm).

3.4 Pollen analysis

The CONISS cluster analysis of the complete pollen assemblage resulted in a separation into four phases, each representing a different regional vegetation composition (Fig. 5).

3.4.1 Phase I-P (2528–1816 cal yr BP)

At the beginning of Phase I-P (340–292.5 cm), forests were dominated by Fagus sylvatica (pollen: 24 %–44.5 %) (Fig. 5). P. sylvestris, Betula undiff., Alnus undiff., Abies alba, and P. abies also constituted a high proportion within the arboreal pollen. The latter two species significantly increased their percentage approaching 1810 cal yr BP. This time interval was characterized by high fire activity events as indicated by high CHAR_micro (max 6035–38 139 $\mathrm{particles}{\mathrm{cm}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) values and the presence of Neurospora, as well as Gelasinospora ascospores between 2500–2300 cal yr BP (Stivrins et al., 2019). Towards the end of the phase, Sphagnum began to increase, indicating a shift to moss-dominated peat.

3.4.2 Phase II-P (1816–1092 cal yr BP)

During the second phase, the forests were dominated by F. sylvatica (18 %–54.5 %) (Fig. 5). A. alba and P. abies were the main components of coniferous forests. A significant decline in F. sylvatica occurred between 1280–1210 cal yr BP. During this phase, A. alba and Quercus reached their highest proportion in the forest. Based on indicator pollen counts, human impact in this phase was the weakest along the entire paleorecord. At the end of this phase, crop introduction in the region was observed, mirrored by an increase in Cerealia pollen share. The stable conditions prevailed at that time in the peatland, which was dominated by Sphagnum and C. vulgaris as well as other Ericaceae species.

3.4.3 Phase III-P (1092–366 cal yr BP)

Throughout this zone, arboreal pollen declined from 97.5 % to 77.5 % between 1090–570 cal yr BP (Fig. 5). A substantial decline in late successional species such as F. sylvatica and C. betulus was observed, especially at the end of this phase. Along with the decrease in the proportion of these species, pioneer trees such as P. sylvestris, Betula, and Corylus avellana constituted an increasing proportion among the woodlands. During this phase, a constant share of cultivated indicators (mostly Cerealia undiff. and Secale cereale) was also recorded, especially from 740 cal yr BP. At the same time, a sharp increase in CHAR_micro, as well as coprophilous fungi taxa, was observed. This corresponded with a sharp decrease in the proportion of Sphagnum (740–570 cal yr BP), the values of which increased again at the end of this phase. At the same time, the disappearance of K. deusta was noticeable during this phase.

3.4.4 Phase IV-P (366 cal yr BP–present)

The major deciduous trees that previously formed the stand gradually withdrew from the site, as manifested by decreasing percentages of F. sylvatica (from 18 % to 4.4 %), Quercus, and Corylus avellana (Fig. 5). P. abies and P. sylvestris reached their highest proportions in the forest (13.5 %–58 % and 21 %–32.5 %, respectively), whereas A. alba seemed to decline completely, as evidenced by the disappearance of its percentage share.

3.5 Biomarker analysis

The TLE following Soxhlet extraction ranged between 20–66 mg g⁻¹, with an average of 35mg g⁻¹. n-Alkanes were measured with chain lengths from C₁₉ to C₃₃, n-alkanols from C₁₆ to C₂₈, and n-fatty acids from C₁₄ to C₃₂. The total abundance of n-alkanes ranged from 80–1112 µg g⁻¹, n-alkanols from 47–469 µg g⁻¹, and n-fatty acids from 260–4224 µg g⁻¹. Overall, the signatures are indicative of higher-plant source material (Eglinton and Hamilton, 1967). The CONISS cluster analysis identified four phases with boundaries at 35 cal yr BP (12 cm), 809 cal yr BP (105 cm), and 1657 cal yr BP (270 cm) (Fig. 6). These four phases are described in the following sub-sections.

Figure 6 Absolute concentration (µg g⁻¹) values for the most abundant homologues of the n-alkanes (ALK; C₂₉, C₃₁, and C₃₃), n-alkanols (ALC; C₂₄, C₂₆, and C₂₈), and n-fatty acids (FA; C₂₂, C₂₄, and C₂₆) as well as the C₂₃ and C₂₅ n-alkanes and C₁₆ n-fatty acid. On the right is the resulting dendrogram of the CONISS analysis of the complete biomarker values with phase boundaries indicated at 35 cal yr BP (12 cm), 809 cal yr BP (105 cm), and 1657 cal yr BP (270 cm).

Figure 7 Carbon preference index (CPI) and average chain length (ACL) values of n-alkanes (ALK), n-alkanols (ALC), and n-fatty acids (FA) and n-alkane, n-alkanol, and n-fatty-acid ratios. The four phases from the cluster analysis are indicated at 35 cal yr BP (12 cm), 809 cal yr BP (105 cm), and 1657 cal yr BP (270 cm).

3.5.1 Phase I-B (2528–1657 cal yr BP)

At the very beginning of the first phase, the TLE was relatively high but began to decrease around 2272 cal yr BP (320 cm). The most abundant homologues for the n-alkanes and n-fatty acids were consistently C₃₁ (33–508 µg g⁻¹) and C₂₄ (293–1376 µg g⁻¹), respectively (Fig. 6). For the n-alkanols, the most abundant homologues varied between C₂₂, C₂₄, C₂₆, and C₂₈ (14–139 µg g⁻¹). The CPI of the n-alkanes (CPI_ALK) ranged from 6.7–14.6, averaging 12.3 (Fig. 7); that of the n-alkanols (CPI_ALC) from 8.7–12.8, averaging 10.7; and that of the n-fatty acids (CPI_FA) from 3.9–8.2, averaging 6.7. The ACL of the n-alkanes (ACL_ALK) varied from 29.0–31.2, averaging 30.3; that of the n-alkanols (ACL_ALC) from 23.8–25.3, averaging 24.8; and that of the n-fatty acids (ACL_FA) from 24.4–25.1, averaging 24.9 (Fig. 7). Of the n-alkane ratios, P_aq ranged from 0.05–0.33, averaging 0.13; P_wax from 0.71–0.95, averaging 0.88; ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ from 0.02–0.18, averaging 0.07; ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$ from 0.08–0.44, averaging 0.20; ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$ from 0.13–0.59, averaging 0.27; and ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$ from 0.44–1.05, averaging 0.67. While ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$ and P_wax generally decreased throughout the phase, P_aq, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ generally increased (Fig. 7). The n-alkanol ratios ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ ranged from 0.55–1.00, averaging 0.75; $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ from 0.62–2.03, averaging 1.03; and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ from 0.34–0.70, averaging 0.52. Generally, $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ increased throughout the phase, while ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ decreased. The n-fatty-acid ratios ${\mathrm{C}}_{\mathrm{15}}/({\mathrm{C}}_{\mathrm{15}}+{\mathrm{C}}_{\mathrm{16}})$ ranged from 0.04–0.08, averaging 0.05, and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ from 0.72–0.90, averaging 0.81. ${\mathrm{C}}_{\mathrm{15}}/({\mathrm{C}}_{\mathrm{15}}+{\mathrm{C}}_{\mathrm{16}})$ stayed nearly the same through the phase, while ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ increased.

3.5.2 Phase II-B (1657–809 cal yr BP)

In the second phase, the TLE decreased further until it began to increase around 1209 cal yr BP (172 cm). The most abundant homologues for the n-alkanes and n-fatty acids consistently remained C₃₁ (19–317 µg g⁻¹), with one exception of C₂₅ (26 µg g⁻¹) at 1223 cal yr BP (175 cm) and C₂₄ (123–747 µg g⁻¹), respectively (Fig. 6). For the n-alkanols, the most abundant homologues varied between C₂₄, C₂₆, and C₂₈ (6–77 µg g⁻¹). All of the CPI values remained well over 1 (Fig. 7). The average of ACL_ALK decreased slightly to 28.9, while the other ACL average values remained nearly the same. The averages of P_aq, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ increased to 0.35, 0.40, 0.55, and 0.19, while the averages of P_wax and ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$ decreased to 0.70 and 0.53. The averages of the n-alkanol ratios ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ decreased to 0.49 and 0.92, while that of ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ increased to 0.61. The averages of the n-fatty-acid ratios stayed nearly the same. During this phase, many of the proxies followed a similar curve or its inverse, reaching either a peak or a dip near the middle of the phase. These include ACL_ALK, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$, P_aq, P_wax, n-alkanol ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$, and n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ (Fig. 7).

3.5.3 Phase III-B (809–35 cal yr BP)

In the third phase, the TLE increased, though it had a slight dip from 466 to 259 cal yr BP (65–40 cm). Following the dip, the TLE increased again to the same level as at the beginning of the phase. The most abundant homologues for the n-alkanes and n-fatty acids were consistently C₃₁ (55–423 µg g⁻¹) and C₂₄ (136–619 µg g⁻¹), with one exception of C₂₆ (242 µg g⁻¹ at 178 cal yr BP; 30 cm), respectively (Fig. 6). For the n-alkanols, the most abundant homologues varied between C₂₄, C₂₆, and C₂₈ (3–103 µg g⁻¹). All of the CPI values again remained over 1 (Fig. 7). The averages of the ACL values all increased slightly (ALK: 30.3; ALC and FA: 25.1). The averages of P_aq, ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ decreased to 0.15, 0.52, 0.19, 0.31, and 0.07, while the average of P_wax increased to 0.86. The average of the n-alkanol ratio ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ increased to 0.62, while that of $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ and ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ decreased to 0.86 and 0.50. The average of the n-fatty-acid ratio ${\mathrm{C}}_{\mathrm{15}}/({\mathrm{C}}_{\mathrm{15}}+{\mathrm{C}}_{\mathrm{16}})$ increased due to a high peak at 250 cal yr BP (40 cm), then decreased to its previous level and stayed almost the same, and that of ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ decreased to 0.78. During this phase, many of the proxies reached a local maximum or minimum around 345 cal yr BP (50 cm) (Figs. 6 and 7).

3.5.4 Phase IV-B (35 cal yr BP–present)

In the fourth, most recent phase, there are only four samples, so patterns are difficult to identify. The TLE initially dipped from its values in Phase III-B but was relatively high in the uppermost sample (43 mg g⁻¹). The most abundant homologues for the n-alkanes and n-fatty acids were consistently C₃₁ (42–75 µg g⁻¹) and C₂₄ (76–263 µg g⁻¹), respectively (Fig. 6). For the n-alkanols, the most abundant homologues varied between C₂₂, C₂₄, C₂₆, and C₂₈ (8–43 µg g⁻¹). The CPI values all remained above 1 (Fig. 7). The average ACL values remained almost the same. The average values of the n-alkane ratios stayed nearly the same, except for ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$, which increased from 0.52 to 0.69 (Fig. 7). The averages of the n-alkanol ratios ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and $({\mathrm{C}}_{\mathrm{22}}+{\mathrm{C}}_{\mathrm{24}})/({\mathrm{C}}_{\mathrm{26}}+{\mathrm{C}}_{\mathrm{28}})$ increased to 0.76 and 0.96, while the averages of the other n-alkanol and n-fatty-acid ratios stayed almost the same.

Table 1 Phases identified with CONISS analysis of different proxies and major changes between each phase.

4 Discussion

4.1 Temporal differences in trends across molecular proxies

The CONISS analysis identified different phases across all of the proxy records, though there were similar trends across most of the records (Table 1). The elemental data were split into five phases: 2528–2113, 2113–1569, 1569–1151, 1151–809, and 809 cal yr BP–present (Fig. 3). The C and N concentrations, $\mathrm{C}/\mathrm{N}$ ratio, and δ¹³C values are reflective of both the source material and the environmental conditions, as well as post-deposition decomposition. Therefore, although there are similarities between the patterns of the elemental data and the other proxy records, they do not match exactly.

The macrofossil data were clustered into four phases: 2528–2251, 2251–1671, 1671–64, and 64 cal yr BP–present (Fig. 4). Compared with the elemental data, the CONISS analysis identified one fewer significant zone, and while the first two phase boundaries of each proxy are within 140 cal yr BP of each other, the others do not align.

The pollen data were clustered into four phases: 2528–1816, 1816–1092, 1092–366, and 366 cal yr BP–present (Fig. 5). The pollen data reflect both a local and a regional vegetation signal, so while there are similarities between the phases based on the pollen and those based on the elemental and macrofossil data, they do not line up exactly.

The biomarker data were also clustered into four phases: 2528–1657, 1657–809, 809–35, and 35 cal yr BP–present (Fig. 6). The first three phases in the biomarker data are the most similar to those of the pollen data, while the final phase is the most closely aligned with that of the macrofossil data. The biomarker data reflect the source material and post-deposition decomposition, so some offsets to the other proxy records are expected.

Because macrofossils are conventionally seen as the most reliable proxy reflecting in situ vegetation and peat development (Birks and Birks, 2000, e.g.,), we used the phases based on the macrofossils to describe the peatland development.

4.2 Development of the Beerberg peatland

In Phase I-M, 2528–2251 cal yr BP, the vegetation is characterized by a dominant presence of E. vaginatum (Fig. 4), along with a relatively low $\mathrm{C}/\mathrm{N}$ ratio and very negative δ¹³C values (Fig. 3), which is typical for fen peat or transitional peat underlying bog peat (Kuhry et al., 1992; Jones et al., 2010). High fire activity during this phase was indicated by the maximum count of macrocharcoal pieces (charred wood, moss stems, Calluna leaves) (Fig. 4) occurring in this interval as well as increased CHAR_micro (Fig. 5) and the presence of carbonicolus ascospores, Neurospora and Gelasinospora (Shumilovskikh and van Geel, 2020). It is likely that some of the fire events were a result of anthropogenic activity, indicated by the presence of Plantago lanceolata or ruderal communities like Artemisia or Chenopodiaceae during this phase (Fig. 5). Fires in the bog area could have initiated the peat development as initial (fen) peat could form on wet ground following fire, developing into transitional peat and then bog peat, as has been seen at other sites (e.g., Tuittila et al., 2007; Gałka et al., 2019). The pollen analysis also indicated that Fagus sylvatica was the dominant arboreal species during this time period (Fig. 5), as was also determined by Lange (1967). Although F. sylvatica is considered a fire-sensitive species (Tinner et al., 2000), it records maximum occurrence in this phase soon after the fire. This could also signify that the fire negatively affected the presence of A. alba and P. abies in contrast to F. sylvatica. The disappearance of Neurospora and Gelasinospora combined with the rapid decline in CHAR_micro corresponds with the development of the Sphagnum community on peatland. This may also indicate wetter conditions on peatlands, especially during 2500–2230 cal yr BP. The increasing $\mathrm{C}/\mathrm{N}$ ratio and δ¹³C values towards the end of the phase (Fig. 3), as well as the initial incidences of Sphagnum (Figs. 4 and 5), indicates a shift towards more ombrotrophic conditions (Wang et al., 2015). The estimated accumulation rate during this phase (Fig. 2) is also relatively low compared with the rest of the core, likely due to the lower proportion of mosses (Stivrins et al., 2017) as it is the paludification phase, the initiation of the peatland. The biomarker data also support this interpretation, with P_wax, P_aq, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ all indicating an initial high proportion of longer-chain-length n-alkanes usually deriving from vascular plants (Fig. 7). Through the phase, the relative abundance of the shorter C₂₃ and C₂₅ n-alkanes (Fig. 6) increases, indicating the beginning of the peat development and an increase in the proportion of Sphagnum mosses (Baas et al., 2000; Pancost et al., 2002; Bingham et al., 2010). Of the other ratios, the n-alkanol ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ generally decreased, while the n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ increased, potentially indicating an increase in n-alkanol C₂₄ and n-fatty acid C₂₄, which have previously been measured as dominant chain lengths in some Sphagnum species (e.g., Ficken et al., 1998a; Corrigan et al., 1976). Additionally, the TLE and individual biomarker absolute concentrations were relatively high compared with later phases in the core (Fig. 6). This aligns with previous research finding that, compared with vascular plant species, Sphagnum mosses produce lower outputs of lipids (Pancost et al., 2002); therefore, as the proportion of Sphagnum increases, the TLE and biomarker concentrations decrease.

Following the shift to ombrotrophic conditions, in the second phase, 2251–1671 cal yr BP, S. fuscum was dominant and the main peat-forming plant (Fig. 4). There was also a fairly steady presence of Ericaceae rootlets, though no evidence of E. vaginatum. At the beginning of this phase, there were also remnants of S. medium/divinum species as well as Polytrichum. In the pollen data, the percentage of Sphagnum spores increased rapidly, along with that of P. abies (Fig. 5). During this phase, the $\mathrm{C}/\mathrm{N}$ ratio and δ¹³C values began to increase (Fig. 3), indicating wetter conditions during rapid peat growth and low decomposition (Loisel et al., 2010; Kuhry and Vitt, 1996). The accumulation rate thus also began to increase during this phase (Fig. 2). In the biomarker data, the n-alkane ratios, indicating an increase in the proportion of shorter chain lengths (C₂₃ and C₂₅), continued to increase due to the rise in Sphagnum. This was also the case for n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ (Fig. 7). The TLE and absolute concentrations of individual biomarkers decreased from the first phase (Fig. 6), as expected with the increase in Sphagnum.

In the third phase, 1671–64 cal yr BP, S. fuscum was still dominant in the macrofossils, but a steady presence of E. vaginatum returned along with S. medium/divinum and Polytrichum in the later half of the phase (Fig. 4). The $\mathrm{C}/\mathrm{N}$ ratio and δ¹³C values reached their peak during this phase, indicating wet conditions and low decomposition (Loisel et al., 2010; Kuhry and Vitt, 1996). Within this interval, Phase II-P (1816–1092 cal yr BP) of the pollen assemblage is contained (Fig. 5). During this period, the presence of the parasitic fungus Kretzschmaria deusta – an indicator of tree fungal infections, often found on deciduous trees, especially on Fagus sylvatica – was recorded (Wilkins, 1934). Among other mountain sites, its presence was mostly associated with a higher proportion of F. sylvatica and other broad-leaved trees (Czerwiński et al., 2020; Kołaczek et al., 2020). The increased presence of this fungus in the past may have been related to stronger herbivore presence and/or coppicing practices and/or grazing damages (Latałowa et al., 2013; Karpińska-Kołaczek et al., 2014). However, the exact interactions between fungal infection and other disturbance factors in the past are not fully understood (Kołaczek et al., 2020). In the case of the Beerberg site, K. deusta was recorded during the period of the lowest human impact and seemed to correspond with the development of forests with a higher role of Carpinus betulus. During this period, CHAR_micro values were very low, which suggests a decline in fire activity or even a lack of fires near the study site. The $\mathrm{C}/\mathrm{N}$ ratio and δ¹³C values generally also began to decline (Fig. 3), which in an ombrotrophic peat indicates drier conditions (Loisel et al., 2010) and increased decomposition (Kuhry and Vitt, 1996) and could be related to human impact, such as drainage of the peatland. In the forest, P. abies expanded (especially from 520 cal yr BP) and possibly replaced F. sylvatica in cleared areas (Fig. 5). This pattern has previously been observed at other sites in Germany, such as the Black Forest (Rösch, 2000; Gałka et al., 2022b). Deforestation was related to human impact, as evidenced by an increase in indicators typical of open landscapes, including meadows and pastures, such as Poaceae and Plantago lanceolata, as well as ruderal habitats (mostly Rumex acetosa/acetosella type, Artemisia, Chenopodiaceae, and Brassicaceae undiff.). Further evidence of human impact results from the steady proportion of cultivated species, as well as the increase in CHAR_micro and coprophilous fungi taxa. During this phase, the biomarker records showed more variation than the macrofossils. The TLE and absolute concentrations were relatively low through the phase until around 809 cal yr BP when they began to increase; there was another dip in the TLE and most individual concentrations (excluding the n-alkanes C₂₃ and C₂₅) from around 466 to 259 cal yr BP, increasing again in the final part of the phase (Fig. 6). The fluctuations in TLE are likely related to a decrease and subsequent increase in Sphagnum, as well as an increase in the proportion of Ericaceae species and other vascular plants (Pancost et al., 2002) in the latter half of the phase. The biomarker ratios follow a similar pattern, with many reaching either a maximum or a minimum point around 1200 cal yr BP, including P_aq, ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{23}}+{\mathrm{C}}_{\mathrm{29}})$, ${\mathrm{C}}_{\mathrm{25}}/({\mathrm{C}}_{\mathrm{25}}+{\mathrm{C}}_{\mathrm{29}})$, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$, as well as the n-alkanol ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ (Fig. 7). These ratios also showed a local maximum or minimum around 345 cal yr BP. This is also seen in the δ¹³C curve as well as in a color change in the peat core itself (Fig. 3). All of these proxies appear to be indicating that there is first a decrease in mosses and water levels beginning around 809 cal yr BP, followed by an increase from around 500 to 345 cal yr BP, and then a subsequent decrease. The earlier decrease around 809 cal yr BP, implying lower water table levels and drier conditions, could be related to the warmer Medieval Climate Anomaly (MCA; ca. 900–1400 CE) (Luterbacher et al., 2016). The following increase in water table levels or surface moisture around 500 cal yr BP could be caused by colder, wetter conditions as a result of the Little Ice Age (LIA; ca. 1300–1850 CE), as similar changes during this interval have been noted in other peatland records from Germany (e.g., Barber et al., 2004) as well as Poland (Marcisz et al., 2020). Following this brief increase, P_aq decreased again, while P_wax increased until around 35 cal yr BP (Fig. 7). This is likely linked to the increased human activity and drainage that occurred at the peatland in the 19th and 20th centuries. As previously mentioned, peatlands throughout Europe have exhibited drier conditions in the most recent times (Swindles et al., 2019).

In the fourth phase, 64 cal yr BP–present, Sphagnum medium/divinum was dominant, and there was an increased presence of E. vaginatum (Fig. 4). This could indicate that the peatland is no longer pristine as the shift in dominant Sphagnum species could be related to pollution through increasing dust deposition (Gałka et al., 2019, 2022a, b) and further dry conditions, related to drainage as well as the warming climate (Swindles et al., 2019). However, in the pollen analysis, the decrease in human impact indicators following 35 cal yr BP could be evidence of the conservation practices introduced on or near the peatland (Fig. 5). This phase coincided with again more negative δ¹³C (Fig. 3) due to the peat here being in early stages of decay (Loisel et al., 2010) and elevated N concentration, potentially caused by the increase in more fen-like vegetation including E. vaginatum (Kuhry et al., 1992). However, the high N concentration could also result from high atmospheric N deposition throughout the last few decades (Ackerman et al., 2019). The biomarker signature in this phase seemed to differ primarily based on a higher abundance of the C₁₆ homologue of the n-fatty acids (Fig. 6). As this corresponds to the current acrotelm of the peatland with active peat formation and high decomposition, the higher abundance of microbial-derived biomarkers is a result of higher microbial activity (Ficken et al., 1998b).

Figure 8 Selected proxy curves from the macrofossil record (horizontal bars: S. fuscum, S. medium/divinum, E. vaginatum), pollen record (filled-in lines: Sphagnum, F. sylvatica, P. abies), and the biomarker and elemental records (lines; n-alkane P_aq, n-alcohol ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$, n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$, and δ¹³C). A composite photo of the core is also shown with the phases from the CONISS analysis of the macrofossils briefly described. The estimated time frames of the MCA and LIA are highlighted on the biomarker and elemental curves in red and blue, respectively.

4.3 Added value of biomarker analysis

Both the biomarker and the macrofossil analyses are reflective of the local vegetation within the peatland, while the pollen analysis reflects both the local and the regional vegetation. Timing differences between local and regional vegetation shifts have previously been reported in studies using biomarkers and pollen analyses (e.g., Jansen et al., 2013).

Despite potential discrepancies in the cluster-analysis-derived phases, the biomarker results illustrate a similar story of peatland development. The primary driver behind the phase changes in the biomarker data appeared to be the abundance of C₂₅ n-alkane. As C₂₅, along with C₂₃, is known to be highly abundant in Sphagnum species, the phases resulting from the CONISS analysis most likely represent changes in the input of Sphagnum to the peat (Fig. 6). Comparing the biomarker results with those of the macrofossils and pollen, the increase in Sphagnum did generally correlate well (Figs. 4, 5, 6, 8). Moreover, the increase in Sphagnum and related biomarkers was paralleled by higher $\mathrm{C}/\mathrm{N}$ ratios and enriched δ¹³C values (Fig. 3), indicative of more ombrotrophic conditions (Wang et al., 2015). The n-alkanol ratio ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and the n-fatty-acid ratio ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ also roughly followed these trends. While they have not previously been linked to Sphagnum mosses as indicative ratios, future peat biomarker studies could investigate this further.

One notable difference shown in the biomarker measurements as opposed to the macrofossils was the behavior of P_aq, P_wax, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$, as well as the n-alkanol ratio ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and the n-fatty-acid ratio ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$ in Phase III-B (Fig. 7). The generally decreasing trend in P_aq and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ at the beginning of this phase implies lower water table levels and drier conditions, which could be related to the warmer Medieval Climate Anomaly (MCA; ca. 900–1400 CE) (Luterbacher et al., 2016). Then from around 500 to 345 cal yr BP, P_aq began to increase, while P_wax decreased. This shift in conditions is also reflected in the δ¹³C values (Fig. 8) and in an abrupt increase in Sphagnum in the pollen data but is not reflected in a notable change in the macrofossils.

Additionally, all of the ACL measures reach a minimum at this same interval (Fig. 6), indicating a higher input of the shorter-chain-length homologues, which would also correspond with less input from vascular plants. Furthermore, the ${\mathrm{C}}_{\mathrm{23}}/{\mathrm{C}}_{\mathrm{25}}$ ratio reached a minimum in this interval, indicating a potential shift in Sphagnum species (McClymont et al., 2008). Following this brief increase, the P_aq decreased again, while the P_wax increased until Phase IV-B or 35 cal yr BP (Fig. 7). This is likely linked to the increased human activity and drainage that occurred at the peatland in the 19th and 20th centuries (Swindles et al., 2019). The indication of P_aq, P_wax, and ${\mathrm{C}}_{\mathrm{23}}/({\mathrm{C}}_{\mathrm{27}}+{\mathrm{C}}_{\mathrm{31}})$ of a change in local conditions that are not reflected in the macrofossils shows that biomarkers can be a valuable additional proxy to use in paleoenvironmental studies.

4.4 Comparison with other local and regional records

The developmental trajectory of the Beerberg peatland was similar to other peat bogs in Germany and other parts of central and northern Europe. The fen-to-bog transition occurred around 2250 cal yr BP, which has been recognized as the beginning of a shift to cooler and/or wetter periods in previous studies in Ireland (Barber et al., 2003), England (Wimble, 1986), and Denmark (Aaby, 1976). The next major transition seen in the core record was the Medieval Climate Anomaly (MCA) and Little Ice Age (LIA), the exact time of which varies across European paleoclimate records (Wanner et al., 2022). Here at the Beerberg peatland, we found evidence of warming in the biomarker measurements and the pollen record during the period ca. 1150–1300 CE, as shown by the decrease in Sphagnum in the pollen record as well as corresponding changes in the biomarker ratios. This was followed by a period of likely cooling from ca. 1450–1700 CE, shown by the return of Sphagnum in the pollen record and corresponding fluctuations in the biomarker ratios. We believe that these periods correspond to the MCA and LIA. This roughly overlaps with results from previous studies: for example, Mauquoy et al. (2002) identified the LIA as lasting from ca. 1450–1600 CE at Lille Vildmose in Denmark and ca. 1460–1600 CE at Walton Moss in England, Barber et al. (2004) identified the LIA as starting between ca. 1250–1350 CE at Dosenmoor in Germany and Svanemose in Denmark, and De Vleeschouwer et al. (2009) found evidence of the LIA from ca. 1200–1800 CE. Following the LIA, the Beerberg pollen record showed that the forests of the region transitioned from primarily F. sylvatica forests to P. abies forests, with an increase in P. sylvestris in the last few centuries as well. This pattern has been seen in many other northern and central European forests (e.g., Svobodová et al., 2002). The final vegetation transition at Beerberg beginning in the late 19th century, with an increase in Ericaceae species and E. vaginatum as well as a shift in dominant Sphagnum species, is likely related at least partially to human influence such as drainage of the peatlands and change in forest composition in the catchment. This has also been noted in other European records (Barber et al., 2004). Generally, the development of the Beerberg peatland confirmed patterns seen at other European sites, though locally transitions occurred with some offsets because of the exposed location of the sequence at the mountain ridge. As such offsets have been well documented, this shows that studying local archives are the most informative of local climate and environmental shifts. To the best of our knowledge, this study is the first of its kind covering the last 2 millennia in a peatland in central Germany, and its location in Thuringia enables future application to research identifying climate boundary shifts in central Europe (Breitenbach et al., 2019).

5 Conclusions

We found that the peatland itself, with Sphagnum fuscum as the dominant peat-forming species, did not undergo much vegetational change following its initial development. This stability persisted even amidst notable shifts in forest composition, from being beech- to spruce-dominated, and increased anthropogenic land use within the region. In the last couple of centuries, the pristine plant population of the peatland was disturbed, most likely by dust deposition and hydrological changes, as we were able to glean from the elemental and biomarker analyses despite the relative homogeneity of the macrofossil analysis.

This study has further demonstrated opportunities for biomarker analysis to contribute meaningfully to paleoenvironmental investigations. Specifically, the biomarker record serves as an independent confirmation of the trends found in the pollen and macrofossils, providing more confidence in the vegetation reconstruction, especially on a species level. The n-alkane ratios provided more precise information about fluctuations in the local conditions of the peat bog, pointing to potential influences from regional climate shifts that underpin the observed changes in vegetation from the pollen and macrofossil data. We also identified two more potential diagnostic biomarker ratios for Sphagnum abundance: the n-alkanol ${\mathrm{C}}_{\mathrm{22}}/{\mathrm{C}}_{\mathrm{24}}$ and n-fatty acid ${\mathrm{C}}_{\mathrm{24}}/{\mathrm{C}}_{\mathrm{28}}$. Additionally, the fact that multiple parameters, such as various sources of organic matter and processes like degradation and preservation of organic matter, can be assessed confirms the high potential for biomarker applications in peat records. However, the numerous molecular proxies derived from biomarker composition are often difficult to interpret independently, requiring certain expertise and the assessment of biomarkers from several compound classes to gain supportive data for certain interpretations.

Consequently, to increase the effectiveness and efficiency of biomarker analyses, a more systematic approach is required, aiming at integrating biomarker compound classes more holistically. The vast majority of biomarker indices and ratios have been developed for n-alkanes, neglecting other compounds such as n-fatty acids and n-alkanols. While the n-alkanes were the most useful compound class for our study, it would be beneficial if other compounds and potential diagnostic ratios were investigated and applied more systematically, as we have done here.

Appendix A: Supplementary data

Radiocarbon dates

Table A1 Radiocarbon dates of selected organic remains. The sample from a depth of 7.5 cm was excluded from the age–depth model.