Landscape reconstructions for Europe during the late Last Glacial (60–20 ka BP): a pollen-based REVEALS approach

Vegetation change during the Last Glacial period in Europe plays a crucial role in better understanding the ecosystem dynamics response to abrupt climate change. Yet, most quantitative reconstructions of land-cover primarily focus on the Holocene and aim to disentangle the impact of anthropogenic and climatic stress on the vegetation. Here, we present temporally continuous land-cover reconstructions from Europe for the latter half of the Last Glacial period (60–20 ka BP) using the “Regional Estimates of Vegetation Abundance from Large Sites” (REVEALS) model. The pollen-based REVEALS model uniquely factors in plant-specific parameters, such as relative pollen productivity (RPP) and pollen fall speed to model pollen dispersal and thus provides more accurate representation of past vegetation cover than simple relations of fossil pollen counts. We compiled a total of 61 datasets from Europe and its bordering regions to model land-cover estimates across 60 time steps in 1000 year increments (Kern, 2024b). By grouping the 38 analysed taxa into 5 land-cover types (LCTs), we simplify the interpretation of our results, demonstrated in three selected periods during the Last Glacial: Greenland Interstadial 14 (GI-14), Greenland Stadial 9 (GS-9), and the Last Glacial Maximum (LGM). These periods provide insight into stadial-interstadial vegetation variability as well as extreme glacial conditions, which seem to play a fundamental role in the demographic developments of Palaeolithic hunter-gatherers. Additionally, we compare the REVEALS land-cover estimates to raw pollen data, provide REVEALS standard errors, and discuss the reliability of our results as well as potential avenues to further improve the reliability of REVEALS estimates. To facilitate the use and interpretation of our data for a wide scientific audience, we developed the browser-based application PALVEG (https://oakern.shinyapps.io/PALVEG/, last access: 22 October 2025), which requires no prior programming experience and dynamically generates maps based on user input.

1 Introduction

The quantification of terrestrial vegetation dynamics in response to climatic change remains poorly understood as it necessitates the conversion of fossil pollen data to tangible reconstructions of the past vegetation cover (Dallmeyer et al., 2023). Yet, terrestrial vegetation is one of the major components of the Earth's terrestrial biosphere and a major driver in the global carbon cycle (Masson-Delmotte et al., 2021). Through biogeochemical feedback processes, their role as a mitigation agent is of particular importance in the context of climate change (Anderegg et al., 2020; Harper et al., 2018; Williamson, 2016). Making sense of the past and developing future adaptation and mitigation strategies greatly depend on our understanding of the interaction between climate and vegetation – an integral part of Earth system models (ESMs) employed to generate climate projections (Fisher et al., 2018; McDowell et al., 2020). ESMs often rely on dynamic global vegetation models (DGVMs) that predict changes in land-cover, productivity, biomass, and/or carbon storage as a function of static bioclimatic variables that are provided by climate models (e.g., Shao et al., 2018; Lu et al., 2018; Hopcroft et al., 2017). Evaluating and improving DGVMs requires appropriate data.

An ongoing endeavour is the comparison of DGVM outputs with palaeoproxy records of the past vegetation composition (Shao et al., 2018; Woillez et al., 2013). The primary goals are to evaluate the output DGVMs and, in a next step, to inform DGVMs and thereby improve the accuracy of their predictions. Pollen records from peat deposits, lake sediments, marine sediments, and to a lesser extent from outcrops and archaeological layers represent the most promising empirical data of past vegetation cover (Dallmeyer et al., 2023). Pollen records provide direct evidence of the vegetational composition and are widespread across terrestrial landscapes. However, the pollen-vegetation relationship is non-linear and often biased by taxon-specific characteristics that influence pollen productivity, dispersal, and deposition (Andersen, 1970; Prentice, 1985). Moreover, boundary conditions, such as CO₂-levels, atmospheric patterns, and soil properties vary through time and may influence these taxon-specific parameters (Ladeau and Clark, 2006; Wayne et al., 2002) thus hindering a direct comparison with DGVM outputs.

Despite these proxy-specific drawbacks, quantitative vegetation reconstructions have a long history in the geosciences (Chevalier et al., 2020). Over time, different methods have emerged and improved upon, such as biomisation (e.g., Prentice et al., 1996; Binney et al., 2017; Li et al., 2025), the Pseudobiomisation Method (PBM; Fyfe et al., 2010; Woodbridge et al., 2014), modern analogue techniques (MATs; e.g., Zanon et al., 2018), and landscape reconstructions (REVEALS; Sugita, 2007; Githumbi et al., 2022; Serge et al., 2023; Schild et al., 2025). The different approaches have been discussed, compared, and evaluated previously (Chevalier et al., 2020). Here, we focus on the Regional Estimates of Vegetation Abundance from Large Sites (REVEALS; Sugita, 2007). Among the aforementioned techniques, REVEALS uniquely factors in pollen productivity and dispersal properties, yielding land-cover percentages for each individual taxon. Syntheses of regional to continental REVEALS-specific parameters are continually expanded and improved upon (e.g., Githumbi et al., 2022; Mazier et al., 2012; Schild et al., 2025; Serge et al., 2023; Wieczorek and Herzschuh, 2020). In the recent years, REVEALS land-cover simulations were generated on global (Schild et al., 2025) and continental scales for Eurasia, China, and North America (e.g., Roberts et al., 2018; Githumbi et al., 2022; Li et al., 2023). However, previous works have mainly focused on the response of terrestrial ecosystems to anthropogenic stress (e.g., agriculture and deforestation) during the Holocene (11 700 yr BP to present).

The Holocene is a period of relative climate stability and low-amplitude climate fluctuations in comparison with the preceding glacial period (Mayewski et al., 2004). In order to focus on the impact of climatic stressors on the regional vegetation, we direct our attention on the Last Glacial period (ca. 115 000–11 700 yr BP), which features marked climatic fluctuations (Fig. 1) with limited human impact on the environment (Brovkin et al., 2021; Fletcher et al., 2010). The Last Glacial encompasses a wide range of climate states: Following the relatively stable interglacial climates of early Marine Isotope Stage 5 (MIS 5e), phases of climate instability (remainder of MIS 5a–d and MIS 3) alternated with relatively stable glacial climates (MIS 4 and MIS 2; Fig. 1). Particularly the climate in the Northern Hemisphere during MIS 3 is characterised by a marked climate instability that manifested in several periods of abrupt warming that permeated the otherwise glacial climatic conditions (e.g., Rasmussen et al., 2014; Corrick et al., 2020). These periods were termed Greenland Interstadials (GIs) after their discovery in Greenland ice cores (Svensson et al., 2006), while periods in-between GIs are termed Greenland Stadials (GSs), accordingly. While these events are relatively well-studied in the marine realm, their impact on terrestrial palaeoenvironments is still largely unknown (Fletcher et al., 2010; Moreno et al., 2014), especially in Central Europe (Britzius et al., 2024; Kern et al., 2022; Woillard, 1978), a fact that is commonly ascribed to the scarcity of available data. Although the climatic boundary conditions today are different, the response of terrestrial ecosystems to abrupt warming events is of particular interest in the context of modern climate change (Masson-Delmotte et al., 2021). The vegetation response to GI-GS variability raises several questions regarding adaptation, resilience, and extinction/repopulation rates (Alley et al., 2003; O'Neill et al., 2017), the possibility of local/regional plant refugia (Tzedakis et al., 2013; Willis et al., 2000; Willis and van Andel, 2004), and also seems to have high relevance for the large-scale demographic developments of Palaeolithic hunter-gatherers (Maier et al., 2022).

Figure 1 Top panel: map of Europe depicting the location of all pollen archives included in this study. The vertical and horizontal grey lines subdivide the study area into four quadrants, which are referenced in the discussion. See Table A1 for a detailed legend and references. Last Glacial Maximum ice-sheet extent redrawn after Ehlers et al. (2011). Bottom panel: Greenland oxygen isotopes (Rasmussen et al., 2014) for the last ca. 125 000 years. Red intervals are referred to in the results and discussion section (Sect. 3). Approximated presence of different human species (Staubwasser et al., 2018) is indicated by the coloured bars. Marine Isotope Stages (MIS) after (Lisiecki and Raymo, 2005).

2 Methods

2.1 Pollen data – sources and preparation

A total of 61 pollen records with data for the 60–20 kyr BP interval were compiled from study sites located in Europe and bordering regions (Fig. 1). The majority of those datasets were obtained from the NEOTOMA database (Williams et al., 2018; Kern, 2025), the ACER pollen and charcoal database (Sánchez Goñi et al., 2017), and PANGAEA (de Beaulieu, 2010; Britzius et al., 2024; Kern, 2021; Koutsodendris et al., 2023a; Miebach et al., 2019a; Sadori et al., 2018). Additional datasets were gathered directly from the respective publications and by personal communication with other researchers. A complete table listing all sites and references can be found in Table A1.

Most datasets were either provided with a robust age-depth model in the corresponding publication or database. Records without a robust chronology were tentatively correlated with well-established records based on the stratigraphic placement of these records in a broader context. Datasets from marine records were excluded due to their overrepresentation of Pinus pollen grains (Heusser and Balsam, 1977). For the purpose of this study, raw pollen counts (compared to e.g., pollen percentages) are required (see Sect. 2.3 for details). However, some datasets only provide pollen percentages and the original data could not be obtained. To maximise the amount of records included in the study, percentage data from these sites were recalculated to a baseline counting sum of 300, which is considered the minimum for statistically robust results and thus often used as a target during palynological analysis (Weng et al., 2006).

Data preparation included the removal of empty and duplicate entries and excluding samples outside of the desired age interval to speed up the REVEALS calculation. Subsequently, all data were aggregated into 1000-year time bins following Githumbi et al. (2022) and Trondman et al. (2016). Since not all pollen records cover the entire target interval, the availability of suitable datasets for each time bin varies with time (see Fig. B1). Data aggregation reduced the input dataset from 5231 to 1632 samples. Afterwards, the harmonisation of nomenclature and pollen morphological types for all datasets was required (Birks et al., 2023). Taxa harmonisation was performed in R using and adapting the taxon harmonisation table supplied by Githumbi et al. (2022).

2.2 REVEALS

The REVEALS model is based on the r-value model (Davis, 1963) and subsequent iterations (Andersen, 1970; Parsons and Prentice, 1981) and its development is described in detail in earlier studies (e.g., Sugita, 2007). In essence, REVEALS reconstructs land-cover estimates for all taxa, factoring in pollen counts, basin parameters, such as basin type (i.e., lake or bog) and basin radius (assuming a circular basin), and plant physiological parameters (pollen productivity and pollen fall speed). Early versions of REVEALS were limited to large lakes (>50 ha) based on the premise that sediments from large lakes accurately represent regional vegetation assemblages (Sugita, 2007). Subsequent development stages of REVEALS then successfully demonstrated that REVEALS estimates from small lakes (<50 ha) and bogs (both small and large) are comparable to those of a large lake in the same region (Trondman et al., 2016). The substantially higher standard deviation for these data can be offset by including multiple small sites, which is consequently highly encouraged (Trondman et al., 2016).

2.3 REVEALS input parameters

The REVEALS model estimates land-cover percentages using raw pollen counts from palaeoclimate archives, supplied on either a depth or an age scale. In addition to pollen data, the REVEALS model requires additional taxa-specific parameters, such as the relative pollen productivity (RPP) estimates and the respective standard deviation as well as the fall speed of pollen (FSP). RPP and FSP values have been previously compiled on continental to global scales (Githumbi et al., 2022; Serge et al., 2023; Wieczorek and Herzschuh, 2020). Moreover, REVEALS land-cover estimates using variable RPP values have been previously validated and optimised through comparison with satellite-based land-cover reconstructions (Hansen et al., 2013; Schild et al., 2025). Here, we use combined RPP and FSP values from Githumbi et al. (2022) and Serge et al. (2023) that are optimised for Europe, but briefly also discuss differences in land-cover estimates by using others sets of RPP and FSP (Table 1). For maximum comparability, we use the same taxa across for all model runs and fill missing values as needed from Githumbi et al. (2022) and Serge et al. (2023).

Table 1 Fall speed (FSP), relative pollen productivity (RPP) including standard deviation (SD) and associated land-cover type (LCT) for all taxa included in this study. We list three sets of RPP values: a synthesis of RPP values for Europe, modified from Githumbi et al. (2022) and Serge et al. (2023); a global compilation (listed are values for Europe, or, if not available, the Northern Hemisphere; Wieczorek and Herzschuh, 2020, dataset v2), and optimised RPP values using a satellite-based validation and optimization routine (Schild et al., 2025). All values for RPP and SD (RPP) are relative to Poaceae (bold).

The REVEALS model has previously been used to estimate anthropogenic land-cover change during the Holocene and includes taxa of significance for this specific research question, such as Cerealia (cereals) and Secale (rye), which are indicators of agricultural practices. Since agriculture only emerged during the early to mid-Holocene in Europe (Price, 2000; Zeder, 2011), we grouped these taxa together with other Poaceae (grasses) in our analysis. Plantago lanceolata (plantain), while also a common indicator for agricultural practices, is grouped together with other taxa from the genus Plantago, as pollen from this family naturally occur in pollen records from Europe throughout the Last Glacial period. The respective RPP value is calculated as the mean of the available values for different Plantago species (Githumbi et al., 2022).

Pollen data from different archives are typically heterogeneous in terms of nomenclature and the level at which taxa are identified. Tree taxa are most commonly identified on the genus (e.g., Picea) or species level (e.g., Picea abies). However, herbs and grasses are typically identified on the family level (e.g., Ericaceae and Poaceae), owing to the sheer diversity of these groups (the Poaceae family contains more than 10 000 species in over 700 genera), rendering the exact identification using microscopic analysis almost impossible. To include as many data as possible, it was necessary to harmonise the pollen data. Given the limited number of datasets to begin with, we resorted to simplifying all datasets to the most common denominator (family level) rather than excluding datasets to overcome the heterogeneity in pollen data. For some taxa, RPP values for boreal/temperate and Mediterranean specimens exist (Githumbi et al., 2022). Since we focus on the Last Glacial, i.e., a period when mostly cold and dry climates prevailed in Europe, we chose the boreal/temperate values, as we expect Mediterranean taxa to play a subordinate role under glacial conditions compared to their boreal/temperate counterparts. As a result, our final RPP and FSP table consists of entries for 38 taxa (Table 1).

Additional parameters used by REVEALS are basin radius (in m), basin type (i.e., lake or bog), and the geographical location of the archive (in latitudinal and longitudinal coordinates). This information was either taken directly from the respective publications, included in database metadata, or roughly estimated using satellite imagery. During the Last Glacial, peat bogs were much less common than during the Holocene and many lakes transformed into bogs as the climatic conditions became warmer and wetter (e.g., Füramoos, Germany). Hence, we set the basin type to “lake” for all sites with no information, unless the respective publications indicated that a bog was present throughout most of the history of the archive (e.g., Tenaghi Philippon, Greece).

REVEALS also requires two independent parameters: wind speed and the maximum extent of the regional vegetation (Z_max). Following a previous evaluation of REVEALS estimates combined with empirical data, wind speed and Z_max were set to 3 m s⁻¹ and 100 km, respectively (Trondman et al., 2015).

2.4 Implementation of REVEALS

The REVEALS model was implemented using the R package LRA (Abraham et al., 2014). The REVEALS function within LRA distinguishes between data from lakes and bogs and accounts for the taxon-specific FSP and RPP (including standard deviation). LRA offers two different dispersal models: a Gaussian plume model (used in this study; Sutton, 1953) and a Lagrangian stochastic model (Theuerkauf et al., 2016).

Previous studies have mostly implemented a grid-based approach to combine records across 1°×1° grid cells in Europe for a combined mean cover estimate and standard error for all taxa. Owing to the scarcity of datasets available in Europe during the pre-Holocene, such an approach is not feasible, as the typical distance between two pollen archives exceeds 1° in either direction. Hence, we chose not to apply this technique to our dataset and provide site-specific estimates, rather than regional (grid-cell wide) estimates. To provide maximum clarity, we specify which sites are considered “small” (i.e., <50 ha) and advocate for a careful interpretation of REVEALS estimates from these sites.

The REVEALS output generates land-cover estimates (including standard errors) for all taxa from all samples provided. For the purpose of this study, we provide REVEALS estimates from all sites in 1000-year increments, with a ±500-year data aggregation window to account for age-depth model uncertainties. As a result, some sites may not be present during all time slices, depending on their temporal resolution and/or potential gaps in the respective datasets (see Figs. 3–5).

2.5 Land cover types

In the biomisation approach (e.g., Prentice et al., 1996; Binney et al., 2017; Bigelow et al., 2003) each taxon is assigned a plant functional type (PFTs), based on their ecological and climatic preferences. Subsequently, groups of PFTs form different biomes, such as desert, steppe, various types of tundra and forests. To simplify interpretation, PFTs can also be grouped into “Super-PFTs” (Binney et al., 2017) as an alternative to biomes. The REVEALS model only includes a limited number of taxa due to the limited availability of reliable RPP estimates. Therefore, we adapted a simplified approach introduced by Githumbi et al. (2022) and group taxa into five land-cover types (LCTs): Open land (which includes all types of grasses and herbaceous taxa), conifers, cold deciduous trees, temperate deciduous trees, and Mediterranean trees. Additionally, we provide estimates for total trees, calculated as 100 % – Open land. See Table C1 for a detailed list of all taxa included in each LCT. The choice of LCTs is based on the vegetational composition during the glacial conditions of the Late Pleistocene and aims to simplify and illustrate broader changes in the vegetation rather than highlight a single taxon.

2.6 Mapping

The results of the REVEALS land-cover reconstructions are presented as maps in time slices of 1000-year steps for all taxa and all LCTs. Figures 3–5 illustrate exemplary maps during full glacial (LGM, 23 ka BP), interstadial (GI-14, 52 ka BP), and stadial (GS-9, 39 ka BP) conditions for the open land LCT. We selected these time slices for their significance in palaeoclimate research. All maps share common features to facilitate their interpretation: Land-cover percentages (REVEALS-based) and pollen percentages are displayed on a gradual colour scale with varying limits to accommodate the vastly different percentages between LCTs and taxa; the differences between REVEALS estimates and pollen counts are also given in percentages on a diverging colour scale (strong colours denote a strong difference, faded colours denote a small difference, white signals that both values are roughly equal); only sites with data for a given time interval (selected time ± 500 years) are displayed; icon size denotes whether a site is classified as small (<50 ha) or large (>50 ha). Sites marked in red denote unreliable data (see Sect. 4.1 for more details).

2.7 R-based application PALVEG

To facilitate the access to our REVEALS land-cover reconstructions of Europe during the Late Pleistocene, we provide an easy-to-use interactive map-generating application: PALVEG, available at https://oakern.shinyapps.io/PALVEG/ (last access: 22 October 2025). PALVEG is a browser-based application using the shiny R package (Chang et al., 2023) and thus does not require any prior coding experience. PALVEG offers easy access to all REVEALS land-cover reconstructions included in this study in the form of dynamically generated maps based on the inputs chosen. PALVEG will be continuously updated with new REVEALS land-cover estimates as more pollen datasets and REVEALS parameters (e.g., RPP and FSP for additional taxa) become available.

PALVEG allows users to define a time slice between 75 and 15 ka BP, in intervals of 1000 years, and chose a taxon or LCT they are interested in. PALVEG will then visualise the associated REVEALS land-cover percentages of all available datasets for the defined time on a map of Europe (10° E–60° W; 30–70° N). PALVEG can also display the pollen data as raw count-based percentages for the same interval. To highlight the differences between both methods, users also have the option to display the deviation of the REVEALS model from raw pollen counts: Blue (positive) values indicate that the REVEALS model attributes a higher land-cover estimate than the count-based estimate. Accordingly, red (negative) values suggest that the REVEALS model assumes a lower land-cover percentage compared to the pollen data. Lastly, PALVEG offers to display both REVEALS land-cover estimates and pollen data aggregated into 1000-year time bins (default) as well as data based on the sample closest to the chosen time slice (within ±500 years).

We acknowledge the need for intercomparability between different REVEALS land-cover estimate maps, which necessitates a uniform colour scale (i.e., ranging from 0 %–100 %). Yet, land-cover percentages for many taxa are well below e.g., 10 % across Europe and individual differences are masked by a uniform scale. For maximum clarity, PALVEG therefore offers two modes to display REVEALS land-cover estimates: A fixed scale (0 %–100 %) and a dynamic scale that changes depending on the maximum land-cover percentages for a given taxon or LCT. This allows for maximum comparability between individual outputs, while also offering a more detailed and readable display option.

3 Results and discussion

In the following, we discuss the spatio-temporal evolution of land-cover occupied by trees in Europe during the Last Glacial using the RPP synthesis described in Sect. 2.3 (Table 1). For simplicity, we refer to total tree populations, which are predominantly characterised by conifers, while other LCTs play a minor role. In general, we observe a sharp, Europe-wide decrease in land-cover of trees starting with the end of MIS 5 (ca. 70 ka BP; Fig. 2a). During the subsequent MIS 4 (70–60 ka BP), a minimum in tree cover is observed in all regions with the exception of NE Europe, where the land-cover of total trees remains at a moderate level between 30 %–40 %, indicating an open woodland vegetation. During early MIS 3 (ca. 60–50 ka BP), tree populations recover, particularly in NW and SE Europe, in response to the climatic amelioration during prolonged Greenland Interstadials (e.g., GI-16, GI-14, and GI-12). In contrast, tree populations decline in NE Europe or are unaffected (SE Europe) in response to the climatic change. The opposite pattern emerges during the later part of MIS 3 (ca. 40–30 ka BP), when a recovery of tree populations can be observed in NE Europe accompanied by a synchronous, but transient, increase of tree populations in SW Europe. MIS 2 (28–15 ka BP) is generally defined by a subdued variability compared to MIS 3. Land-cover percentages of tree populations in NW and SE Europe are close to their respective minimum for the region, whereas NE and SW Europe signal slightly elevated levels of tree cover.

Figure 2 (a) Evolution of the mean REVEALS land-cover percentages for all LCTs through time (75–15 ka BP). The data have been categorised in four quadrants split by latitude (>46 and <46° N, ca. North and South of the Alps, respectively) and longitude (>10.5 and <10.5°E, roughly separating the Western Mediterranean from the Central and Eastern Mediterranean). (b) Site availability by quadrant for each time step.

Figure 3 Results for the LCT “Open land” during GI-14 (52 ka BP). (a) REVEALS land-cover percentages. (b) Combined open land pollen percentages (see Table 1). (c) Difference between REVEALS land-cover percentages and the combined pollen data. Colours indicate a higher (blue) or lower (red) land-cover representation compared to pollen data. (d) Same as (a), but data considered unreliable are highlighted in red (see Sect. 4.1 for details).

Overall, the regional vegetation across Europe is characterised by a division into two major regimes: in NW and SE Europe, the vegetation pattern follows that of other proxy data such as NGRIP δ¹⁸O (Fig. 1; Rasmussen et al., 2014): Periods of warmer and more humid conditions coincide with elevated levels of land-cover occupied by trees and vice versa. This pattern is related to large-scale shifts in atmospheric circulation patterns and thus moisture availability (Florineth and Schlüchter, 2000; Ludwig et al., 2016). Similar observations in speleothems and loess records have been ascribed to shifts in the North Atlantic Oscillation and associated storm tracks responsible for moisture distribution (Luetscher et al., 2015) and dust transport (Schaffernicht et al., 2020) across Europe.

For NE Europe, the opposite pattern emerges. It appears that during the coldest and driest intervals (MIS 4 and MIS 2) of the Last Glacial, tree populations remain at a moderate level and instead decline during the warmer and more humid period of MIS 3. However, we ascribe this pattern to data scarcity rather than climatic variability. Particularly in NE Europe, very few datasets are available (Fig. 2b) and they span a large latitudinal range (from 46 to 67° N). During mid-MIS 3, only limited data are available and thus the mean forest cover is strongly influenced by records from the high latitudes (e.g., Lake Yamozero, site 2, Fig. 1) or in close vicinity to ice sheets (e.g., Nesseltalgraben, site 21, Fig. 1), which tend to signal a very open landscape. During subsequent intervals (i.e., late MIS 3), more records from, e.g., the Carpathian Mountains, are included and thus shift the mean towards higher forest covers. Note that despite the overall site availability in NW Europe being similar or even lower compared to NE Europe (Fig. 2b), these sites span a much smaller geographical area and climatic range. As a consequence, the presence or absence of individual sites has a reduced impact on the regional means.

In the following, we focus on three critical intervals of the Last Glacial, using the LCT “open land” as proxy for vegetation openness: Stadial-interstadial variability, represented by GS-9 (39 ka BP) and GI-14 (52 ka BP), respectively, and the Last Glacial Maximum (LGM, 23 ka BP). To facilitate the comparison, we treat these climatic events as singular time slices centred around their respective mid-points (given in parentheses). Then, we compare the output of the REVEALS model with pollen percentages and highlight some of the major differences in the perception of palaeoenvironments. Lastly, we discuss the reliability of REVEALS estimates in the context of the REVEALS model design and its fundamental ideas.

3.1 Stadial-interstadial variability

For the comparison of interstadial to stadial land cover, we chose the longest interstadial (GI-14) and stadial (GS-9) to fully capture the differences in land-cover composition. Both climatic events are well-documented in numerous records and chronologically well-constrained. Also, shifts between stadial and interstadial conditions and their associated ecological gradients may have had considerable impact on the large-scale distribution of regional populations of hunter-gatherers during the Upper Palaeolithic (Maier et al., 2024).

In contrast to the interstadial conditions of GI-14, average open land percentages increased to 80 % during Greenland Stadial 9 (GS-9, 39 ka BP) and are highest along the northern fringes of the Alpine glaciers (Fig. 4). This increase is primarily due to a rise of Poaceae to 32 %, while Cyperaceae (15 %) and Asteraceae (16 %) remained on similar levels as compared to GI-14. Interestingly, Poaceae percentages not only increased during GS-9, but the geographical distribution of Poaceae spread to most regions of Europe. Notable exceptions are in the Iberian Peninsula, the Alpine Foreland, and the Carpathian Mountains, where a strong local signal of Cyperaceae indicated the presence of wetlands. Asteraceae and Apiaceae are primarily found in the Central Mediterranean and South-West Asia, locally reaching more than >40 % of combined land cover.

Figure 4 Results for the LCT “Open land” during GS-9 (39 ka BP). (a) REVEALS land-cover percentages. (b) Combined open land pollen percentages (see Table 1). (c) Difference between REVEALS land-cover percentages and the combined pollen data. Colours indicate a higher (blue) or lower (red) land-cover representation compared to pollen data. (d) Same as (a), but data considered unreliable are highlighted in red (see Sect. 4.1 for details).

Overall, tree populations in Europe collapsed during GS-9. While the western Balkans and Eifel region show a declining presence of tress, locally constrained areas in the Carpathian Mountains (87 % trees) and the Baltics (up to 49 % trees) appear to have served as refugia for Picea abies, Betula, and Pinus during glacial periods (Fig. 4). However, caution is advised for the Carpathian data (Fig. 4d). Small, isolated communities of temperate deciduous trees are present in Northern Africa (13 %), the South-West Asia (13 %) and the Eastern Mediterranean (6 %) and comprise the remainder of land covered by trees during GS-9.

Across all sites, the REVEALS-derived mean land-cover estimate of “open land” during Greenland Interstadial 14 (GI-14, 52 ka BP) was 74 % (Fig. 3). The primary constituents of open land vegetation during GI-14 show distinct geographical patterns: Poaceae (26 %) and Cyperaceae (17 %) are the dominant taxa in Central Europe and the western Mediterranean, whereas Asteraceae (15 %) is found in South-West Asia and the Mediterranean coastal regions. Among other taxa, isolated populations of Apiaceae and Chenopodiaceae are present in South-West Asia and Ranunculaceae reach >10 % in Central European wetlands environments. Other herbs only play subordinate roles.

Although an open steppe-tundra vegetation dominated the European landscapes during GI-14, two distinct forested areas can be identified: The western Balkans and to a lesser extent the Eifel and eastern Baltic regions (Fig. 3a), where conifers constitute the majority of the vegetation. However, caution is advised regarding the reliability of some of these data (Fig. 3d). In contrast, the Italian Peninsula shows a mixed temperate woodland, consisting of deciduous Quercus, Fagus sylvatica, Acer, and Abies alba at the time. North of the Alps, forests composed of boreal tree taxa such as Pinus and Betula constitute a substantial portion of trees while other tree taxa are virtually absent.

Stadial-interstadial variability is primarily characterised by an increase in open-land percentages during stadials and a decrease in open-land percentages during interstadials, but probably also by a more pronounced differentiation of the onset of the vegetation period in the southern parts of Europe during stadials (Maier et al., 2024). Such increases in open-land percentages have led to a southward displacement of the tree line in western Central Europe, while limited tree population in eastern Central Europe persist. During interstadials, decreases in open-land percentages are accompanied with a northward expansion of the tree line. These observations are in line with the interpretation of proxy records across Europe (Fletcher et al., 2010; Landais et al., 2022, Tzedakis et al., 2013). However, the scarcity of available pollen records for GI-9 from Northern and North-eastern Europe substantially hinders our ability to further inquire. Additional pollen records are required to investigate the spatio-temporal framework of tree-line recession in a north-easterly direction during interstadials in more detail.

The arrival and spread of Palaeolithic hunter-gatherers between 43 to 40 ka BP in Europe (Shao et al., 2024) occurred during a phase of pronounced interstadials (ca. GI-12 to GI-9) and therefore relatively mild and humid climate conditions (Rasmussen et al., 2014). This climate state is also reflected in increased tree cover percentages, particularly in Eastern Europe (i.e., the Balkans and along the Danube), where climatically suitable east-west corridors may have temporarily opened up for human dispersal across Europe (Shao et al. 2024). Here, light woodland or a mosaic of open and forested vegetation was prevalent particularly in SE Europe (Fig. 2a). Such ecotones provide ideal conditions for a high faunal biodiversity and are thus attractive habitats for hunter-gatherers. During stadials, the tree line shifted southwards and the landscape in the higher latitudes opened up. At the same time, phenological gradients related to the greening of the landscape in spring became more pronounced in the southern parts of Europe. The gradual cooling of interstadials likely led to a shift in largescale phenological patterns, presumably affecting the spatial distribution of populations and providing incentives to move, for instance, into the more southern areas of the Iberian Peninsula (Maier et al., 2024), or foster retreat to local niches (Timmermann, 2020). Combined with decreasing temperatures and moisture availability, much of Central and Eastern Europe might have become unfavourable for hunter-gatherers during pronounced stadials, such as Heinrich Event 4, and populations may have decreased (Shao et al. 2024). There are hints that Neanderthals had smaller mobility ranges than anatomically modern humans and occupied more fragmented habitats (Timmermann, 2020). Together, this may have fostered sensitivity to climate change during MIS 3 and its implications for ecosystems (Yaworsky et al., 2024), contributing in the long run to their disappearance.

3.2 Last glacial maximum

REVEALS land-cover estimates for the LGM are characterised by uniformly high percentages of open land vegetation (88 %) across most of Europe. The decrease in forested area by 40 % from GS-9 (or 48 % from GI-14) demonstrates the substantial impact of the increasingly cool and dry climates on the European vegetation going from MIS 3 to MIS 2. Similar to the previous intervals, the main constituents are Poaceae (31 %), Cyperaceae (18 %), and Asteraceae (23 %), although their geographic distribution varies substantially. Poaceae is mainly found in southern Central Europe and the Western Mediterranean region, while Cyperaceae is primarily found at higher latitudes (>46° N) and Asteraceae at lower latitudes (<46° N). Similar to previously discussed intervals, Apiaceae is most prominent in South-West Asia (>30 %). Substantial amounts of Artemisia (up to 37 %) and Chenopodiaceae (up to 10 %) are scattered across Europe.

Although tree populations are close to their minimum extent during the LGM (Fig. 5a), two locally constrained glacial refugia exist, where trees constitute more than 30 % land cover: In the Carpathian Mountains, Picea abies and Pinus demarcate an open coniferous woodland vegetation. Similarly, strands of boreal forests composed of Picea abies, Pinus, and Betula prevail in the Balkans. This presence of substantial tree populations highlights the existence of glacial tree refugia with suitable microclimates during the LGM.

Figure 5 Results for the LCT “Open land” during the LGM (23 ka BP). (a) REVEALS land-cover percentages. (b) Combined open land pollen percentages (see Table 1). (c) Difference between REVEALS land-cover percentages and the combined pollen data. Colours indicate a higher (blue) or lower (red) land-cover representation compared to pollen data. (d) Same as (a), but data considered unreliable are highlighted in red (see Sect. 4.1 for details).

Among other sites, the total tree-cover during the LGM is very low and reaches a maximum of 25 %. Nonetheless, geographical patterns are discernible from the REVEALS estimates. Small populations of Pinus are present in a punctuated, SW–NE belt spanning from Iberia to Russia, while being absent in the Central and Western Mediterranean as well as North-western Europe. Cold deciduous trees and temperate deciduous trees are scattered across Europe but only occur locally in low percentages (typically <10 %). Noticeable amounts of Mediterranean tree taxa during the LGM are only attested in Greece (Carpinus orientalis; 1.5 %).

3.3 Comparison of REVEALS estimates and pollen data

The direct comparison between REVEALS land-cover estimates and raw pollen percentages highlights the differences between both approaches (Figs. 3c, 4c, and 5c). Expectedly, an overrepresentation of trees in pollen records (and consequently much lower REVEALS land-cover estimates) can be observed (Fig. 6). Accordingly, grasses and herbs are underrepresented in pollen data and achieve much higher REVEALS land-cover estimates. However, on a taxon level, this overarching pattern may not always apply and individual taxon-specific land-cover percentages ultimately depend on the overall vegetational composition. REVEALS estimates may be lower or higher for certain taxa than the respective pollen percentages (e.g., for Poaceae and Cyperaceae). In most cases, this is associated with the taxon-specific RPP values in relation to the overall RPP values of all other taxa present at the site (Table 1). Notable outlier to this pattern is Abies alba, which despite relatively high RPP values is considered underrepresented by the REVEALS model, owing to the high fall speed (FSP) of Abies alba pollen grains and therefore its limited spatial dispersal capabilities.

Figure 6 Comparison of REVEALS land-cover percentages for three critical time intervals, using different datasets of relative pollen productivity (RPP) as discussed in Sect. 2.3 (Table 1). The data are compared to the fossil pollen percentages.

As a result of the disparity between REVEALS data and fossil pollen data, the presence of forests and (open) woodlands may be exaggerated when interpreting fossil pollen data. However, the extent of e.g., the over-representation of tree taxa is highly variable through time and land-cover estimates require careful analysis in order to fully utilise their potential. Overall, these observations are in line with the general conception among palynologists and palaeobotanists that fossil pollen data do not reflect the palaeovegetation one-to-one, which ultimately led to the inception of models such as REVEALS (Andersen, 1970; Prentice, 1985; Sugita, 2007).

3.4 The impact of RPP on REVEALS land-cover estimates

In Sect. 2.3 we have discussed the existence of various compilations and syntheses of RPP and FSP values. While our study used a synthesis of RPP data from Githumbi et al. (2022) and Serge et al. (2023), we want to briefly address the differences in REVEALS land-cover estimates using other sets of RPP values (Table 1). All available datasets are established using high-quality estimates of pollen productivity compiled from research of the last decades and all RPP datasets focus on slightly different aspects in terms of regional (e.g., Mediterranean) to global scale, environmental variability, and methodology. Hence, there is no best dataset for all purposes.

In general, all RPP datasets yield a vegetation composition that indicates an overrepresentation of trees in fossil pollen data and consequently much higher land-cover estimates for “open land” (Fig. 6). However, there are several distinct differences between the RPP datasets. First and foremost, the results for the synthesis used in this study and the compilation by Wieczorek and Herzschuh (2020) are strikingly similar, with the latter dataset generating slightly higher land-cover estimates for the different tree LCTs. In contrast, land-cover estimates using the optimised RPP (Schild et al., 2025) result in very low percentages of tree LCTs. In detail, this primarily affects “conifers”, influenced by the vastly higher RPP of Pinus and to a lesser extent Picea abies, which both affect land-cover estimates negatively. Other tree LCTs are in good agreement with those of other RPP datasets. Differences among herbaceous taxa vary across RPP datasets, but are less significant on regional scales.

4 Data reliability

The reliability of REVEALS land-cover estimates is derived from the standard errors (SEs) that are provided alongside each REVEALS model run (Serge et al., 2023; Trondman et al., 2015). The relative size of SEs is considered a measure of the quality of the REVEALS output. REVEALS land-cover SEs mainly rely on the quantity and quality of the available pollen datasets. This includes both the number of available datasets as well as their individual counting sums for each sample and the classification into large (>50 ha) and small (<50 ha) lakes and bogs (Githumbi et al., 2022; Sugita, 2007). Additionally, the standard deviation of the REVEALS parameter RPP is factored into the estimation of the REVEALS SEs.

A key aspect of previous land-cover reconstructions using the REVEALS model for the Holocene period is to present the results in 1°×1° grid cells, which aggregate data from all sites within a grid cell, based on the type and size of a basin (e.g., Githumbi et al., 2022; Serge et al., 2023; Trondman et al., 2016). The idea is that the quality of REVEALS estimates for a grid cell increases with the number of large lakes present within each grid cell (Sugita, 2007). Later studies have demonstrated that pollen data from small sites (<50 ha) will also improve REVEALS estimates, although to a lesser extent (Fyfe et al., 2013; Trondman et al., 2016). Hence, the quality of vegetation estimates has been related to the number and basin size of the archives where pollen records originate from. For the Last Glacial period, this approach is not feasible due to the low geographical cover of suitable (i.e., chronologically well-constrained) pollen archives in Europe (Fig. 1). Based on the current data compilation, only three 1°×1° grid cells contain more than one pollen archive from the Last Glacial period. Therefore, we provide point-based land-cover estimates that represent the local to regional vegetation within a maximum radius of 100 km around the site (as defined in the REVEALS parameter Z_max).

4.1 Estimation of uncertainty

To evaluate whether REVEALS estimates are reliable, we follow the protocol of Githumbi et al. (2022): If the standard deviation of the REVEALS estimate is equal to or larger than the REVEALS estimate itself, then the data is considered unreliable. This excludes entries with zero counts (and thus a REVEALS land-cover of 0 %), which are considered reliable by definition. In total, 39 % of all non-zero REVEALS land-cover estimates are considered unreliable (Table D1). However, a substantial amount of unreliable data is related to taxa that only appear in very low counts (e.g., Urtica, Rubiaceae, Potentilla, and Filipendula), which are expected to have a high standard error (Sugita, 2007) and/or have high RPP standard deviations (e.g., Abies alba and Populus for which no data is considered reliable). In contrast, taxa that are commonly found in glacial climates such as Poaceae (2 %), Cyperaceae (8 %) or Pinus (2 %) are on average considered highly reliable. Establishing LCTs further increases the reliability compared to the individual reliability of most their constituents, e.g., only 2 % of open land and 14 % of total trees datapoints considered unreliable.

Due to the fact that many taxa are absent in a lot of the samples, the perceived reliability of REVEALS land-cover estimates is much higher. For example, only 8 % of all Urtica samples are considered unreliable if all samples are considered. However, 95 % of datapoints are classified as unreliable if zero-count samples are excluded. Hence, we advise caution regarding the use and interpretation of data from relatively uncommon taxa. Instead, LCTs represent a broader spectrum of the vegetation and their data are more reliable and thus less prone to misinterpretation.

4.2 Data aggregation

Previous studies have suggested to aggregate samples from intervals with a high temporal density in order to increase the size of pollen counts and thus minimise REVEALS SEs (Githumbi et al., 2022; Trondman et al., 2016). In this study, samples ± 500 years around the studied time interval were aggregated into a single sample, which was then used to estimate REVEALS land-cover percentages. As a result, the input dataset shrunk from 5231 to 1632 samples, corresponding to an average aggregation ratio of 3.21 to 1 across all samples. However, due to the unevenly temporally spaced nature of palynological datasets, the number of aggregated samples per bin may vary substantially. The number of sites that contain data for each time step are given in Fig. B2.

Overall, data aggregation decreased the number of unreliable REVEALS land-cover estimates from 48 % (raw pollen counts) to 40 % (aggregated pollen counts). The percentage of unreliable data was reduced drastically for almost all taxa and LCTSs. The biggest improvements occur for taxa with low to very low counting sums per sample, such as Carpinus orientalis (39 % more reliable), Buxus sempervirens (20 % more reliable), Acer (18 % more reliable). Exceptions are Poaceae (0.3 % less reliable) and Abies alba (no improvement). Among LCTs, the reliability of estimates for cold deciduous trees, temperate deciduous trees, and Mediterranean trees increases (9 %, 8 %, and 11 %, respectively), whereas the reliability of estimates for conifers decreases slightly (0.3 %). In conclusion, data aggregation has substantially improved to REVEALS land-cover estimates among less common taxa and LCTs, whereas estimates for the dominant taxa and LCTs remain virtually unchanged.

The percentages of unreliable REVEALS land-cover estimates for both aggregated and raw data can be found in Table D1. Based on our findings, we suggest to use the REVEALS land-cover estimates from aggregated data due to the predominantly higher data reliability observed. Nonetheless, to provide maximum clarity, we provide REVEALS land-cover estimates for both the aggregated and raw datasets and also include both datasets in the PALVEG application.

5 Code availability

The REVEALS function is part of the R package LRA (Abraham et al., 2014), available at https://github.com/petrkunes/LRA (last access: 29 April 2024). The full code used in this study is available on Zenodo (https://doi.org/10.5281/zenodo.12625221, Kern, 2024b). It includes data preparation, taxon harmonization, and data analysis.

The application PALVEG is hosted at https://oakern.shinyapps.io/PALVEG/ (last access: 22 October 2025) and the R code is available on Zenodo (https://doi.org/10.5281/zenodo.12624842, Kern, 2024a).

6 Data availability

REVEALS land-cover estimates and REVEALS standard errors generated in this study are available in the PANGAEA database (https://doi.org/10.1594/PANGAEA.973049, Kern et al., 2025). The metadata for pollen datasets are available on Zenodo (https://doi.org/10.5281/zenodo.16812766, Kern, 2025).

7 Conclusions

Over the recent years, the REVEALS model has been successfully implemented in numerous research projects, such as LandClim (e.g., Githumbi et al., 2022) and PAGES LandCover6k (e.g., Harrison et al., 2020), as well as studies (e.g., Roberts et al., 2018; Kaplan et al., 2017). The results have been utilised in a wide array of research disciplines including but not limited to anthropogenic land-cover change, to evaluate archaeology-based reconstruction of land-use in the past (Morrison et al., 2021), to estimate the impact of mining activities (Schauer et al., 2019), or in the context of biodiversity and regional to local scale vegetation dynamics (Marquer et al., 2017).

For the first time, the REVEALS model has been applied to pollen records from Europe during the latest interval of the Last Glacial period, generating continuous (1000-year intervals) land-cover estimates for 38 taxa and 5 LCTs. The Last Glacial is of particular importance for our understanding of abrupt climate change (i.e., GI/GS variability) and the response of ecosystems (e.g., Fletcher et al., 2010; Landais et al., 2022; Rasmussen et al., 2014). Yet, to this day, vegetation dynamics in Europe during the Last Glacial are far from understood and current hypotheses include rapid re-migration of tree taxa in response to climate amelioration, local refugia providing suitable microclimates that sheltered from the harsh environmental conditions, and complex spatio-temporal combinations of the aforementioned mechanisms (Gavin et al., 2014; Hošek et al., 2024; Tzedakis et al., 2013). Moreover, the Last Glacial saw the disappearance of Neanderthals and the dispersal of Homo sapiens into Europe followed by pronounced demographic fluctuations (Hublin, 2015; Mellars, 2004; Schmidt et al., 2021). It has long been hypothesised to what extent climatic fluctuations have impacted these patterns of population dynamics on both temporal and spatial scales (Maier et al., 2022, 2024; Staubwasser et al., 2018; Timmermann and Friedrich, 2016). Here, corridors for the dispersal of Palaeolithic hunter-gatherers into different parts of Europe may have episodically opened up in relation to the opening of the landscape and shifts in ecological gradients. Furthermore, tree populations have persisted in local refugia in Central and Eastern Europe, where favourable environmental conditions could have provided shelter during subsequent stadial climatic downturns. More sophisticated vegetation reconstructions might shine light on open research questions regarding food and prey availability, hunting and subsistence strategies, and demographic developments (Baker et al., 2024; Ordonez and Riede, 2022; Schmidt and Zimmermann, 2019; Seersholm et al., 2020; Vidal-Cordasco et al., 2023).

With the PALVEG application, we aim to facilitate the access to palaeovegetation reconstructions to a wide scientific audience. PALVEG dynamically generates maps of REVEALS land-cover estimate or pollen data for all sites, taxa, LCTs, and time intervals between 75 and 15 ka BP. PALVEG is browser-based, intuitive, and requires no a priori experience in coding or palaeobotany.

8 Outlook

Further improving the REVEALS model highly depends on the availability and quality of plant physiological parameters (RPP and FSP). Researching pollen productivity and dispersal is tedious and relies on numerous field studies and although steady progress has been made over the last couple of years, more RPP and FSP data are required. This is of particular significance for taxa that are presumed to shape the European landscape during glacial times (e.g., differentiating between tree and shrub growth forms for taxa such as Betula and Pinus; Birks, 1968) or those that are well-known for their proportionately high or low pollen productivity (e.g., Salix and Larix; Niemeyer et al., 2015) and are thus misrepresented in pollen records. Hence, the most promising avenue to further improve land-cover estimates is to provide these parameters for more taxa, ideally on the species or genus level, to fully utilise the information that available pollen datasets provide. Specifically for the Late Pleistocene in Europe, more pollen data at a sufficiently high temporal resolution are required to further refine the geographical cover of available data and the inferred palaeoecological reconstructions. Ideally, these new data originate from large lakes (>50 ha) and fill geographical gaps in areas of Europe where the present distance between available archives is particularly large, e.g., Central and Eastern Europe (Fig. 1). A different approach to improve REVEALS-based palaeovegetation models lies optimisation of RPP using remote-sensing techniques, which has already shown promising results and is expect to further improve in the future.

In the face of anthropogenic climate change, generating high-quality palaeovegetation data remains imperative for improving coupled vegetation-climate models. Only then can ESMs provide more accurate predictions of vegetation and ecosystem dynamics in the near future.

Appendix A: Pollen archives used in this study

Table A1 List of all pollen archives, REVEALS land-cover reconstructions and their respective references. “Site ID“ refers to the labels in Fig. 1.

Appendix B: Site availability per time window

Figure B1 Available sites with pollen data for each 1000 year time window. The selection of sites was based on the availability of data between 60 and 20 ka BP (dashed lines).

Appendix C: Taxon classification and land-cover types (LCTs)

Table C1 Classification of plant taxa, genera, and families into land-cover types (LCTs), adapted from Githumbi et al. (2022) for the Last Glacial period.

Appendix D: Data reliability per taxon/LCT before and after data aggregation

Table D1 Comparison of the percentage of unreliable data for all taxa and LCTs before and after data aggregation into 1000 year time bins.