Introduction

The archaeological record can be interrogated on several spatial and temporal scales, ranging from decisions of individuals at a specific time and location up to phenomena that emerge only as the aggregated activities of populations over the course of millennia and across large areas. We consider research at all scales necessary and relevant. Importantly, not all aspects of human behavior are analyzable at all scales. Phenomena such as the timing of system responses to internal and external factors or feedback processes can differ strongly between different scales of observation (Zhang et al., 2004) and it is necessary to adjust the scale of analysis to the scale of the phenomenon of interest (Maier et al., 2022a). The acquisition and redistribution of energy and information, for instance, cannot be organized the same way in small groups and large populations, which is probably why human social systems seem to self-organize in hierarchical, self-similar structures over successive levels of scale (Zhou et al., 2005; Hamilton et al., 2007; Fernández-López de Pablo et al., 2022). While individual intentions and decision making is highly relevant for small-scale processes, they become increasingly less relevant at larger scales (Maier et al., 2022b). In other words, the aggregated actions of tens of thousands of individuals (or more) with at times opposing intensions become visible at larger scales as an emergent property of populations.

This study is concerned with the relation between the spatially and temporally large-scale distribution of populations in Europe during the Upper Paleolithic and how phenological patterns of the vegetation at a corresponding scale might have contributed to this distribution. We have argued earlier that shifts in the timing and intensity of the vegetation period have a certain explanatory power for the large-scale demographic developments during the Upper Paleolithic in Europe (Maier et al., 2022c). Here, instead of looking at dynamics in population sizes, we focus on the distribution of regional populations in Europe and explore their potential relation to gradients in the timing of the vegetation period.

At this large spatial and temporal scale, fine-grained causal chains are difficult to establish and we would already be satisfied if we observe correlations between vegetational patterns on the one hand and the distribution of human populations on the other. However, a tentative causal chain that links patterns in vegetation to the distribution of human populations is as follows: Human populations in the Paleolithic sustained themselves in part from plant food (Piperno et al., 2004; Aranguren et al., 2007; Pryor et al., 2013), which is why there is a direct relation between the seasonal availability of certain plant resources and the presence of hunter-gatherers in the landscape. This relation, however, seems rather weak given that terrestrial animals apparently played the major role in subsistence (Stevens et al., 2010; Bocherens et al., 2014; Drucker et al., 2017). The most important link between the vegetation and human population therefore is indirect and mediated by herbivorous animals moving through the landscape. These movements are partly influenced by phenological gradients of the vegetation. Fresh plants have a particularly high nutritional value and play an important role in the diet of animals after the end of the winter, when the young are born (Ibid.; Esmaeili et al., 2021). Many herbivorous animals thus follow the greening of the landscape in spring, a behavior sometimes called ‘surfing the green wave’ (Merkle et al., 2016; Aikens et al., 2020). The timing and productivity of the vegetation period thus likely had an influence on the presence of animal biomass in the landscape and, in turn, on the human population preying upon it.

Timing and productivity of the vegetation, however, were not stable and likely changed markedly between stadial and interstadial conditions during the Upper Paleolithic (Rasmussen et al., 2014). Our basic assumption for the following reflections therefore is that the timing and productivity of the vegetation period have a certain explanatory potential for the dispersal and distribution of human populations in larger regions and longer periods of time, inasmuch as they provide long-term and large-scale trajectories for shifts in the abundance of vegetational and animal biomass in the landscape. On small temporal scales, these trajectories might be comparatively weak. However, over the course of millennia, also weak differences in incentives for movements might bias the actions of many animal and human individuals towards preferred directions. Down this line, they potentially have an accumulative influence on the distribution of animal biomass and human populations.

We are aware of the fact that differences existed in regional plant composition as well as that movements of herbivores are only partly linked to the greening of the landscape and that differences in diet exists between different species (Guthrie, 1990). Conspecific density of animals or long-term shifts in aridity (e.g., between stadials and interstadials) certainly also were an important factor and clearly differences existed between species (Aikens et al., 2020). We also acknowledge that not all members of a species show the same migratory behavior at all times (e.g., rather stationary reindeer south of the Loire; Fontana, 2017, 2023). Consequently, we do not aim at reconstructing the vegetation period of specific plant communities or seasonal movements of specific animal species in a certain region — not to mention movements of specific human groups or individuals. To avoid misunderstandings, the following points cannot be overstated:

  • We neither claim that the observations discussed below are the only or the most important factor for the large-scale distribution of human populations in the landscape, nor do we intend to downplay the role of social factors.

  • When we talk about vegetation, we neither refer to a specific patch of land with a specific plant community, nor to specific plant species, but to the vegetational biomass as a whole.

  • When we talk about animals, we neither refer to specific animal individuals, nor to animal species, but to the herbivorous biomass as a whole.

  • When we talk about humans, we do not refer to specific groups or individuals, but to the human population in the investigated area.

  • When we show estimates of phenological gradients for vegetation in spring an autumn, the focus is neither on specific seasonal or annual cycles, nor on the exact timing of the onset and end of the growing season, but on large-scale and long-term patterns generally prevailing during spring and autumn over many thousands of years

  • Finally, when we talk about movement, we do not refer to specific routes that animal or human individuals or groups might have taken over a year or a lifetime, but rather to aggregated movements of animal biomass and shifts in the location of human population over millennia at a subcontinental scale.

In this study, we present a case-transferable approach to diachronically estimate the timing of the vegetation period and visualize resulting phenological gradients. We then discuss the implications of our results for two complementary case studies. First, we look at the Aurignacian dispersal of anatomically modern humans (AMHs) in Western and Central Europe. Past research has done much to illuminate the course and timing of the dispersal of AMHs in Europe (Higham et al., 2011; Shao et al., 2021), but little is known about potential factors guiding this movement apart from larger topographic features such as rivers or the sea shore (Conard & Bolus, 2003; Mellars, 2011). Second, we look at the East European Plain (EEP) during the Middle and Late Upper Paleolithic. Here, a large variety of taxonomic units is reported (Noiret, 2009; Sinitsyn, 2010), but little explanation can be found as to why the archaeological record appears so remarkably heterogeneous. With our approach, we shed light on the explanatory potential of large-scale phenological patterns for guiding population dispersal and fostering the formation of distinct expressions of material culture in the absence of topographic divides. We thus focus on one factor among numerous others and potentially more important factors without claiming a decisive or deterministic role.

Materials and Method

Selected Areas of Interest

To explore the explanatory potential of our approach, we selected two case studies that contrast strongly in topography and archaeological background — namely the Aurignacian in Western and Central Europe and the Gravettian and Epigravettian in the East European Plain — and thus allow for a complementary evaluation.

The Aurignacian is connected to the extensive dispersal of anatomically modern humans in Europe roughly between 43 and 33 ka. In comparison to later periods, it displays a rather homogeneous material culture record (Gennai, 2021). The area of Western and Central Europe covers about 3 million km2 and is characterized by a highly structured topography with contoured shorelines and large mountain chains, such as the Pyrenees, Alps, or Carpathians that structure the landscape. The distribution of Aurignacian sites in Central Europe shows a clear spatial structuring, where areas with high site densities, called Core Areas (Schmidt et al., 2021), alternate with areas void of sites (Schmidt & Zimmermann, 2019). The most important Core Areas of the Aurignacian population are located in southwestern France and northern Spain in Western Europe and Lower Austria/Moravia/southern Poland in Central Europe. Smaller populations can be found in the area of the Upper Danube and Belgium.

The East European Plain (EEP) covers about 4 million km2, stretching from the Vistula River in the West to the Urals in the East and from the White Sea in the North to the Black Sea in the South. Contrary to Western and Central Europe, and despite its enormous extent, the topography is rather uniform with comparatively low elevations and numerous rivers and virtually no natural barriers for the movement of animals or humans. The Middle and Late Upper Paleolithic (roughly 33 to 14 ka) in the EEP post-dated the extensive dispersal of modern humans, but a high population dynamic is evidenced by genetic turnovers (Posth et al., 2023). In contrast to the relative homogeneity of the material culture record of the Aurignacian in Western and Central Europe, the Middle and Late Upper Paleolithic in the EEP is characterized by an apparently highly heterogeneous material culture record with several smaller, contemporaneous taxonomic units being reported from different areas, such as the Brynzenian, Prut culture, Lypa culture, and Lower Dniestr culture (Noiret, 2009), as well as the Streletskian, Gorodtsovian, and Spitsynean, the latter two of which are exclusively connected to the locality of Kostenki (Sinitsyn, 2010). The validity of these taxonomic units as representing social phenomena of comparable spatial scale levels and cultural coherence can and must be debated (e.g. Reynolds and Riede, 2019; Riede et al., 2020; Maier et al., 2022b). However, the fact that it is also difficult to subsume the respective assemblages under one or two larger units in a satisfying way indicates that the technological and morphological variability of assemblages in the EEP is comparatively pronounced. Given that interaction usually promotes similarities in material culture items (Boyd & Richerson, 1985), while non-interaction promotes differences (Koerper & Stickel, 1980; Neiman, 1995) and given further the lack of topographic barriers impeding interaction across the EEP, this cultural heterogeneity is particularly intriguing. Within the EEP, we selected 11 areas of interest because of their geographic location and richness or scarcity of Gravettian and Epigravettian sites. At the western fringe of the EEP lies the Polish Jura, a long-term hotspot of Paleolithic research (e.g., Kot et al., 2019). The eastern slopes of the Carpathians and border region between Romania, Ukraine, and Moldova — here referred to as “East Carpathian” — is likewise known for its ample record of important Gravettian and Epigravettian sites (e.g., Noiret, 2009; Steguweit et al., 2009; Anghelinu et al., 2021). The northern Pontic area, in contrast, is characterized by a relative scarcity of sites (Demidenko, 2018), while the regions of Volhynia and the Mid Desna, both located above 50° northern latitude, show an important, but punctuated record (Chabai et al., 2020; Chabai & Vasyliev, 2021). The Kostenki-Borshchevo region provides an exceptionally rich record of stratified sites, some of which have pronounced local and regional idiosyncrasies (Sinitsyn, 2010). Sungir (Bader, 1978), famous for its Gravettian burials (Trinkaus et al., 2014), has been selected because of its northern position. From the northern Caspian areas, no Middle or Late Upper Paleolithic sites have been reported so far. The Urals delimit the investigated areas to the east. In the EEP’s southern part, the presence of cave art has been reported (Dublyansky, et al., 2018). The central and northern parts are of interest, because of their positions at high latitudes, reaching up to the Arctic Circle (Svendsen et al., 2010).

Estimating Parameters of the Growing Season

For our study, we use a script based on Python and R which can be downloaded here: https://github.com/lou-tha/PaleolithicGreening.

The start of the growing season is temperature-driven and, by most definitions, the start of the vegetation period is set to the first of six consecutive days of a mean daily temperature near the surface above 5 °C. Analogously, the end is set to the first of six consecutive days below 5 °C (e.g., Van Meerbeeck et al., 2011; Jiang et al., 2018; cf. Maier et al., 2022c). This is, of course, a general value and some species diverge from it. However, since we do not want to consider specific plants, we will use this established threshold for our estimates. Since mean daily temperatures in the climate models underlying our estimates are averaged over 30 years, increases and decreases are rather steady and temperature reversions are rare. To account for this increased steadiness, we define the beginning of the growing season at the onset of the first three consecutive days above 5°C (i.e., the first three consecutive of the three-day rolling averages that is greater than 5°C). Accordingly, the end is set to the last day before the first three consecutive days below 5°C.

Since temperatures play a vital role in the timing of the vegetation period, estimating diachronic changes requires spatial information on daily temperatures throughout the year. To date, such paleoclimate information for the Weichselian Glacial is available for 21 ka (LGM). The LGM is one of the focus time slices of PMIP3 (Paleoclimate Modelling Intercomparison Project, Phase 3, Braconnot et al., 2012), with a variety of global climate models modeling the climate during the LGM. For this study, the regional climate model WRF (Weather Research and Forecasting, Skamarock et al., 2008) was nested into the coarse gridded MPI-ESM-P (Max Planck Institute for Meteorology Earth System Model in Paleo model, Giorgetta et al., 2013) LGM simulation. Regional paleoclimate modeling allows for an improved representation of the spatial variability of climate variables due to their refined horizontal and vertical grid spacing (Ludwig et al., 2019). The 30-year-long regional climate simulation provides the average daily temperature (2 m above the ground) at 50 km grid spacing for the two case studies (Fig. 1). The WRF model was adjusted to paleoclimatic settings (e.g., ice sheets, lowered sea level, orbital parameters; all based on PMIP3 protocol) for 21 ka to consider the glacial environmental boundary conditions (Ludwig et al., 2017). Atmospheric boundary conditions in the WRF model were updated at 6-hourly intervals based on MPI-ESM-P data. A second adjustment of sea surface temperatures based on the MARGO project (Multiproxy Approach for the Reconstruction of the Glacial Ocean surface, MARGO Project Members, 2009) leads toward a better agreement of the WRF model results with the available proxy data (Ludwig et al., 2017).

Fig. 1
figure 1

Model domains for regional model experiments with the WRF model for the two case studies (a) East European Plain (EEP), (b) Western/Central Europe Figures include LGM topography (in meters above LGM sea level, shaded), additional land areas due to lower sea level and borders of ice sheets considered in the WRF model

For the analysis, we extracted the daily temperature estimates of the WRF model from a raster image between 50 and 65° northern latitude and from 47.5 to 65° eastern longitude for the EEP case (Fig. 1a), and between 35 and 55° northern latitude and from 10° western to 30° eastern longitude for the Western and Central Europe case (Fig. 1b) with a distance between the sampling points of 2.5° in both spatial directions.

In order to expand the estimates of the timing of the growing season to other stadials and interstadials, diachronic temperature records are needed. To date, reconstructed diachronic temperature estimates are available for Greenland and the Black Sea. The Greenland record provides temperature inferred from δ15N of entrapped air of the NGRIP ice core (Kindler et al., 2014). However, it is unclear to what extent these findings are transferable to other areas of Europe (Margari et al., 2020). The amplitude of differences in atmospheric temperature is influenced by the high seasonal fluctuations in Greenland, probably not directly comparable with the more balanced continental situation in the EEP. The Black Sea gravity core (25GC-1) was drilled in 2007 in the south-eastern area of the Black Sea at a depth of 418 m covering the period between 64 and 20 ka with a precision of 240 years (Wegwerth et al., 2015). The sample locations in the Black Sea are thus closer to the EEP study area, but incorporated the buffering effects of water in contrast to differences in terrestrial air temperatures. None of these records are thus optimal for our study. The differences in atmospheric temperatures over Greenland may have been stronger than over most parts of the EEP, while the buffered differences in water column of the Black Sea were certainly smaller than that of the atmosphere over the EEP. Consequently, the differences between interstadials and stadials in the Black Sea-driven scenario must be considered a minimum, while the Greenland-driven scenario might overstate the differences to some degree, thus providing a frame of reference for the EEP. For the Aurignacian, we only use the NGRIP values, because the Black Sea record is unlikely to provide meaningful data for Western Europe.

Both records only provide annual temperature estimates. Differences in temperatures, however, are not evenly distributed throughout the year, but seem to have been more pronounced between December and May and less pronounced between June and November (Flückiger et al., 2008). In order to arrive at monthly differentiated estimates, we take the average diversion of the maximum and minimum modeled monthly differences for Central Europe (Ibid.) to calculate the percentage value of temperature change per month throughout a year (For details see Supplementary 1 and 2). These values serve as factors to distribute the annual differences in temperatures over 12 months for each GI and GS between 25 and 41 ka. The thus-derived values are then added to the temperature estimates for each day in the LGM model (Supplementary Table 1). Because of the monthly resolution, daily temperature at the beginning and end of each month will be slightly over- and under-estimated for the first half of the year and vice versa for the second half, but we consider these divergences negligible given the scale of our study. From the thus-obtained daily temperature values, we estimate the start, length, and end of the growing season (Supplementary Table 2). To ensure an ecological meaningful length, we included all areas with a vegetation period of ≥ 21 days in the following analyses. Based on these estimates, we then determine the strength of phenological gradients between different locations. To this end, we use the length of the offset between the start and end of the vegetation period at two neighboring points in the spatial grid (2.5×2.5°) as a measure of the strength of potential trajectories between these points (Supplementary Table 3). Thus, the larger the offset, the stronger the gradient. The reasoning here is that if greening starts contemporaneously or only with a slight difference between neighboring areas, the difference in the quality of available plant food is similar in the entire area. Therefore, a caloric investment in moving over longer distances would not result in the benefit of gaining food of higher quality. Put another way, if the grass is not greener elsewhere and the quality of plant food in my neighborhood is the same than in 1000 km distance, the pull of the ‘green wave’ is weak and the incentive to move is therefore low, because there will be no gain from going there. It is the time lag in the greening and the resulting difference in food quality that sets the incentive to move over longer distances. If, therefore, greening in one of two neighboring patches is markedly delayed, the pull is stronger and therefore so is the incentive to move. That does not mean that animals would not move in areas with simultaneous greening, but the lack of a clear gradient would likely lead to more random and less directed movements. Since values and gradients are calculated across the entire map section, some nodes are overlying the ocean. These estimates can be used to trace large-scale shifts in the web of gradients, but are, of course, not informative for terrestrial vegetation patterns.

To be sure, our approach to reconstruct daily stadial/interstadial temperatures should not be understood in competition with other model data. Ideally, high-resolution paleoclimatic model data with daily temperature values would be available continuously for all periods and areas. This is, however, unfortunately not the case. We, thus, aimed at a procedure that allows for generating estimates also for periods where no sophisticated climate models are available. If model data is available, the presented method can also be run using the model data directly. Model data always comes with uncertainties. To assess the robustness of our estimate, we compare our results for GI-11 and GS-8, GI-8 and GS-6 to those published by Armstrong et al. (2019). We find that our approach generally returns higher summer and lower winter temperatures. These differences are slight in most parts of the study area and highest in the north-west at the Atlantic Coast (see Supplementary Tables 4 and 5 and Supplementary Figs. 1 and 2). Using the data provided by Armstrong et al. (2019) would thus mainly result in a later onset and later end — and thus a shift — of the growing season. It should be mentioned that although a good fit exists between the chronological sequence of the empirical NGRIP temperature data — on which our estimates are based — and the model data by Armstrong et al. (2019, Fig. 10), there are also clear offsets between the absolute temperature values. In the future, it would be interesting to do an in-depth comparison between different models for the same periods and estimate the uncertainties.

Results

In the following, we present our estimates for vegetational gradients in spring and fall during selected interstadials and stadials. Rather than focusing of the specific aspects, we are interested in the generic stadial and interstadial signals observable in the two case studies. It is important to note that spring and fall estimates are not simply mirroring one another and thus creating cyclical patterns, but result in a net difference of movement incentives. In the long run, these differences in spring- and autumn-time movements will likely lead to a spatial shift of the activity zones of hunter-gathers towards the direction favored by the differences of movement incentives. In the following, we focus on phenological gradients not at annual but centennial or millennial time scales and on the potential incentives for directed movement they might create.

Aurignacian in Western and Central Europe

For the Aurignacian in Western and Central Europe, we present estimates for the interstadials GI-11 (ca. 43 ka) and GI-8 (ca. 37 ka) and the stadials GS-11 (ca. 42 ka) and GS-8 (ca. 36 ka) (Rasmussen et al., 2014), and therefore for the beginning and end of this taxonomic unit. At first inspection of the results (Figs. 2, 3, 4, 5), it stands out that there are more similarities between the two interstadials, and between the two stadials, than there are between periods that are most chronologically proximal to each other (i.e., between GI-11 and GS-11, and between GI-8 and GS-8, respectively). This is not terribly surprising, but it underlines the fact that the generic signal from stadial and interstadial conditions is stronger than individual differences between specific stadials and interstadials. Hence, although some differences can be expected, the observed signals are probably fairly representative also for conditions during the other stadials and interstadials of the Aurignacian.

Fig. 2
figure 2

Western Europe: NGIPR-driven estimates of the growing season and resulting trajectories during spring. Upper: GI 8, lower: GI 11; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period. Core Areas after Schmidt and Zimmermann (2019)

Fig. 3
figure 3

Western Europe: NGIPR-driven estimates of the growing season and resulting trajectories during fall. Upper: GI 8, lower: GI 11; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period. Core Areas after Schmidt and Zimmermann (2019)

Fig. 4
figure 4

Western Europe: NGIPR-driven estimates of the growing season and resulting trajectories during spring. Upper: GS 8, lower: GS 11; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period. Core Areas after Schmidt and Zimmermann (2019)

Fig. 5
figure 5

Western Europe: NGIPR-driven estimates of the growing season and resulting trajectories during fall. Upper: GS 8, lower: GS 11; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull. Triangles show the start of the vegetation period. Core Areas after Schmidt and Zimmermann (2019)

During the interstadial phases (Figs. 2 and 3), strong gradients are observable predominantly along the coastlines and the northern ice margins, while the southernmost parts and the interior of the continent exhibits only weak or no gradients. The areas of weak gradients in the interior of Western and Central Europe correspond to regions in which Core Areas with high site densities are conspicuously absent. There are no gradients into the Iberian Peninsula in spring. To the contrary, there seem to be pull factors from the north-eastern interior to the Pyrenees. In fall, however, and particularly during GI-11, there seem to be some gradients towards the center and the south.

During stadial phases (Figs. 4 and 5), we see a pronounced southward shift of the entire phenological system and strong decrease in the strengths of the gradients. Over the mid-latitudes, the mesh is also much more fragmentary and compartmentalized than during the interstadials. Particularly in spring, there are hardly any gradients forming over Western and Central Europe. This pattern is more pronounced during GS-8 than during GS-11. The gradients in fall generally form a more close-knit mesh. Interestingly, the coastal parts of the Iberian Peninsula now become integrated into the fabric of gradients.

Generally, the dominant direction of the phenological pull in spring is towards the north and in fall towards the south. However, the pronounced orographic differences around the Pyrenees, Alps and Carpathians create countertrends against the general seasonal pull, i.e., southwards during spring and northwards during fall. This phenomenon produces phenological troughs, i.e., sink-areas of phenological pulls that direct the gradients from the surrounding patches towards their center, disconnecting them from the predominant seasonal pull. During the interstadials, they form both in spring and fall over the Franco-Cantabrian region and the Swabian Jura, both areas of high site densities. During the stadials, when the effect seems more pronounced, these troughs form in spring over western France and the Balkans and thus over areas with rather low site densities, but also over the Basque country and in GS-8 over the Czech Republic where site density is high. This can also be observed for the stadial fall situation, when south-western France becomes additionally integrated.

Middle and Late Upper Paleolithic in the East European Plain

Similar to the Aurignacian in Western and Central Europe, for the Middle and Late Upper Paleolithic in the EEP, there are also broad similarities between interstadial patterns on the one hand and stadial patterns on the other, and these generic signals are much stronger than individual differences. For the Middle and Late Upper Paleolithic in the EEP, we have two estimates for each stadial and interstadial based on the temperature record from the NGRIP-core in Greenland and a core from the Black Sea (see Materials and Method section). To exemplify the generic differences as clearly as possible, we selected GI-8 (ca. 37 ka) and GS-6 (ca. 33 ka) (Rasmussen et al., 2014) as examples for interstadial and stadial condition, since of all estimated periods, they show the strongest differences between one another and therefore are most instructive. The maximum differences of temperature in the NGRIP record between GI 8 and GS 6 amounts to more than 20°C, while for the Black Sea record, it is buffered to around 3.5°C. As for the other case study, the results can be considered fairly representative also for other stadials and interstadials. Moreover, since we selected those stadials and interstadials with the largest differences in temperature and because our estimates are temperature-driven, possible deviations will range between the observed values.

The first thing that stands out when looking at the results (Figs. 6, 7, 8, 9) is the marked differences between the estimates for stadials using the temperature proxies from Greenland (Figs. 6 and 7) and the Black Sea (Figs. 8 and 9). The estimates for interstadials are comparable between the Greenland-driven and Black Sea-driven scenarios, and there are also hardly any differences between the interstadial and stadial estimates within the Black Sea-driven scenario. This is particularly the case for the fall estimates. In all these cases, a continuous mesh of gradients is observable for the entire investigated area. With the Greenland-driven scenario, however, these differences are stark and the stadial estimates show a strongly reduced and more fragmented mesh. All estimates for interstadial conditions see the regions of Volhynia, Mid-Desna and Sungir firmly nested within the web of vegetational gradients and even displaying strong pull factors during spring. In the NGRIP-driven estimates for GS-6, in contrast, they are located outside the web of gradients, which are generally rather weak. This observation is even more pronounced for the areas in the Urals as well as the western part of the EEP, in particular for the Polish Jura.

Fig. 6
figure 6

East European Plain: NGRIP-driven estimates of the growing season and resulting trajectories during spring. Upper: GI 8, lower: GS 6; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull. Triangles show the start of the vegetation period

Fig. 7
figure 7

East European Plain: NGRIP-driven estimates of the growing season and resulting trajectories during fall. Upper: GI 8, lower: GS 6; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period

Fig. 8
figure 8

East European Plain: Black Sea-driven estimates of the growing season and resulting trajectories during spring. Upper: GI 8, lower: GS 6; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period

Fig. 9
figure 9

East European Plain: Black Sea-driven estimates of the growing season and resulting trajectories during fall. Upper: GI 8, lower: GS 6; Arrows indicate the direction, color the strength of the gradient (the darker the color, the stronger the pull). Triangles show the start of the vegetation period

Except for the NGRIP-driven estimates for GS-6, the gradients are stronger along the fringes of the EEP, while those in the central area are rather weak, particularly during fall. This finding is consistent with the topographic situation, where the rising elevation of the Carpathian and Ural mountain chains as well as the relative proximity to the periglacial zone of the Fennoscandian ice shield presumably leads to delays in the onset of the vegetation period in comparison to the more central areas of the EEP. A common feature of all estimates is that north-south-directed gradients are generally stronger than east-west-directed ones, reflecting the differences in temperatures related to latitude. In the absence of pronounced orographic differences, and in contrast to our findings from Western and Eastern Europe, there are no countertrends against the general seasonal pull and no phenological troughs.

Roughly around and north of the latitude of the central Urals, the estimated gradients are — except for the periglacial fringe — very weak, or absent in the case of the NGRIP-driven estimates for stadial GS-6. For the spring estimates in the NGRIP-driven interstadial and Black Sea-driven stadial scenario, the latitudinal area of Kostenki-Borshchevo shows a long stretch of weak gradients, separating the latitudes of the northern Caspian and northern Pontic areas from those of the southern Urals.

Discussion

The above-presented results provide new insights into hitherto invisible phenological gradients and boundary situations with explanatory potential for the large-scale distribution and movement of Paleolithic populations. In the following section, we will discuss some aspects relevant to the two case studies against this background.

Aurignacian in Western and Central Europe

There is growing evidence for a pre-Aurignacian presence of anatomically modern humans at different locations in Europe (e.g., Slimak et al., 2022; Hublin et al., 2020), but these people apparently did not contribute genetically to later populations (Posth et al., 2023). The Aurignacian, therefore, is still the earliest period of extensive dispersal of anatomically modern humans throughout Western and Central Europe. In the following, We will discuss how our results conform to different expansion models and hypotheses, such as the coastal route, the Danube Corridor or Ebro Frontier hypothesis.

During Interstadial Conditions, the Coastal Route Provides the Strongest Pull Factors

A prominent hypothesis of the east-west expansion of anatomically modern humans sees two routes as likely corridors for dispersal (Mellars, 2011). On the one hand, the route along the Mediterranean coast, often connected to the Proto-Aurignacian, and, on the other, the route along the Danube (Conard & Bolus, 2003). Our estimates show a continuous chain of fairly strong gradients along the Mediterranean shore, both in spring and fall, forming under interstadial conditions (Figs. 2 and 3). These gradients also continue along the French Atlantic coast all the way up to Britain. This might have explanatory potential for some very early dates from sites with anatomically modern humans in Britain, such as Kents Cavern (Higham et al., 2011). In this case, the coastal dispersal would have continued up to the British Isles (Shao et al., 2021). The central Balkans and central Western Europe, in contrast, show only weak gradients, if at all, thus providing only weak phenological incentives for dispersal. These findings do not support models which see the Danube as a major corridor into Central Europe (Mellars, 2011), but are rather in line with a more critical view (Chu, 2018). It is also in accordance with the observation of a geographic cline in the composition of personal ornament assemblages. The extremes of this cline are located in Austria and Germany, which have no ornament type in common despite their relative geographic proximity and neighboring position along the Danube (Vanhaeren & d’Errico, 2006). During stadial conditions (Figs. 4 and 5), the chain of gradients along the coast disintegrates, as does the mesh over the continent, and phenological troughs become more numerous. Based on these observations, it can be postulated that fast and extensive dispersals may have taken place during interstadials, while stadials may have been periods with a reduced spatial dynamic in population displacement.

During Interstadial Conditions, there were Likely Only Weak Phenological Incentives for an Expansion into the Italian and Iberian Peninsula South of 42°N Latitude

Our results suggest little to no incentive to move into the Italian Peninsula further south than Tuscany during interstadials. Only the northern half is integrated in the mesh of gradients and during spring, pull factors are mainly directed towards the north (Figs. 2 and 3). During stadial phases, however, the system switches. While the northern half becomes rather disconnected, the southern half becomes integrated into the phenological mesh and during spring, Tuscany becomes a phenological trough (Figs. 4 and 5). This might have some interesting implications for the discussion on the distribution and internal dynamics of the Uluzzian in Italy, with its predominantly southern distribution (Peresani et al., 2019).

During interstadials, there also seems to be little phenological incentive for movements into the Iberian Peninsula and an expansion of AMH to the south of the Ebro was probably not fostered. Not only are gradients absent over most of Iberia, but an expansion further south or southwest runs also counter to the northward oriented gradients in spring. During stadials, however, a more or less continuous chain of gradients forms along the shoreline of the entire Peninsula, with particularly strong gradients at the western and southern shore. Towards the center, however, gradients are rather weak. These findings might add a new perspective to the debate about the Ebro Frontier hypothesis (Zilhão, 2000, 2021) and the settlement of Iberia during the Aurignacian.

To avoid misunderstandings, we do not argue that the Italian or Iberian Peninsula were less attractive habitats for hunter-gatherers. To the contrary, longer growing seasons and a probably high plant and animal biomass might have made them quite attractive. The point we wish to make is that, from the perspective of phenological gradients, we see little incentives for movements into these areas during interstadials.

Core Areas are Located at Phenologically Special Positions

Several Core Areas, i.e., areas of high site density, have been identified for the Aurignacian (Schmidt & Zimmermann, 2019). Interestingly, the three most important Core Areas in Franco-Cantabria, southwestern Germany, and Lower Austria/Moravia/southern Poland are located in regions which our analysis identifies as having special phenological settings.

In Franco-Cantabria and southern Germany, the Core Areas coincide with regions of comparatively strong gradients that run counter to the prevailing seasonal movements (Figs. 2, 3, 4, 5). In southern Germany, the structure of pull-factors may have had some isolating effects during interstadial phases, disconnecting this Core Area from the surrounding mesh of gradients, creating a phenological trough. This trough is even more pronounced in stadial phases, when the area is decoupled from the gradient network. In Franco-Cantabria, in contrast, we observe during interstadials a seasonal change between a pulling-apart phenomenon in spring and a pulling-together phenomenon in fall, centered over the Dordogne region. The pulling-apart phenomenon is characterized by gradients pointing away from the center in all directions, while in fall, the pull-gradients are all pointing towards the center. Within such a phenological pattern, stable fission fusion cycles can form that span from the Loire to the Pyrenees and from the Atlantic to the western slopes of the Massif Central. During stadial phases, however, the area becomes largely decoupled from the mesh of gradients, keeping close ties only with the Basque Country, Cantabria and Asturias.

The Core Area in Lower Austria, Moravia and southern Poland is located in an area of weak or no gradients during interstadial periods, but during stadials, it also becomes a central place for pulling-apart and pulling-together phenomena. However, in contrast to Franco-Cantabria, the seasonal order is reversed with pulling-together gradients in spring and pulling-apart gradients in fall. During stadial phases, the region might thus receive influx from populations in the Pannonian Basin.

The East European Plain

Radial Pulls Towards the Fringes of the EEP and Dominant North-South Oriented Gradients Make Networks Forming in Either the Western or the Eastern Part more Likely than Networks Transgressively Spanning the Entire Area

The EEP is characterized by two large drainage systems with a predominantly north-south orientation: the Dnipro and Don Rivers flowing into the Black Sea and the Volga River flowing into the Caspian Sea. Our estimates suggest that incentives for roughly south-north directed movements are generally stronger than those for west-east-directed movements in the entire EEP. That means movements along the river system is fostered more than movement between the river systems, both under stadial and interstadial conditions as well as in both spring and fall. In addition, pull-factors at and towards the fringes of the EEP seem to be stronger than those at or towards the center. This is because the higher altitudes of the surrounding mountains and the periglacial zone cause a delay in the start of the growing season when compared to the central part of the EEP. Although these differences may not be very pronounced on the small scale, over the course of millennia, a slight but constant incentive to move towards the fringes of the EEP might pull populations away from the center. Under these conditions, networks forming in either the western or the eastern part of the EEP seem more likely than networks spanning both parts. Consequently, interaction within the western or eastern part is more likely than between the parts. Given that non-interacting groups are likely to accumulate differences in cultural traits (Neiman 1995), these phenological conditions could have been a factor in the development of idiosyncratic, self-contained manifestations of artifact morphologies that have been archaeologically observed. The low incentives for east-west-directed movements might be part of the explanation why the EEP during the Gravettian and Epigravettian contains so many regional variants without matching equivalents in other parts of the plain, such as the industry of Byki, 18–16 ka (Akhmetgaleva & Demidenko, 2017) at the Seim river in the Dnipro catchment or the Gorodtsovian or Spitsynean in the Kostenki-Borshchevo region at the river Don (Sinitsyn, 2010). This region with its numerous sites from many periods of the Upper Paleolithic is also remarkable in other ways. It is located close to the center of the EEP, where it remains within the mesh of gradients even during pronounced stadial conditions and where, at the same time, pulls towards the west and east are weakest. It may therefore receive population influx during the entire Upper Paleolithic and from both the western and the eastern part, although connections to the west seem to be a bit stronger. This might be part of the explanation why the archaeological record in this region is particularly rich and manifold.

Shifts in Phenological Patterns Between Interstadial and Stadial Periods Might have been a Pace-maker for Human Range Pulsation in the EEP Fostering the Growth and Disintegration of Social Networks

Traditionally, hunter-gatherer presence in the higher mid-latitudes of Europe has been viewed as a more or less continuous phenomenon. Recent research, however, often sees it as subject to range pulsation (Gorodkov, 1986), at least partly influenced by environmental factors (Obreht et al., 2017) and thus as sporadic and punctuated (Maier et al., 2020; Pedersen et al., 2019) and affected by population fluctuations (Maier & Zimmermann, 2017; Maier et al., 2022c).

For the EEP, the non-buffered NGRIP-driven estimates show a strong oscillation in the extent of the mesh of phenological gradients, forming a far-flung web over the entire investigated area during interstadial conditions, while contracting markedly during stadial condition towards the south and east.

During interstadial conditions, all areas of interest discussed in this paper become integrated in the web of gradients, likely supporting network growth and intensification of contact. Also, the area for punctuated presence seems to have expanded markedly to the north. At around 34 – 32 ka, people seem to have been present at the western slopes of the northern Urals close to the polar circle at the site of Byzovaya, 65° northern latitude (Svendsen et al., 2010). The range of available radiocarbon dates cover stadial and interstadial periods (Ibid.). Our estimates suggest that during stadial periods, the vegetation period at these latitudes was likely too short to attract and maintain populations of herbivorous, making a presence during interstadial periods more likely. Also during interstadials, it seems that pulls towards the north start weakening around the latitude of the central Urals.

During stadial conditions, the situation seems to change markedly. Depending on the extent of the cooling, the Urals and areas of Sungir, Mid Desna, and Volhynia are located at the fringes or even outside the mesh of gradients. Phenological incentives for movements between the western part of the EEP and Central Europe seem to be particularly weak, decoupling the Polish Jura from the rest of the network. This suggests that during stadial conditions, phenological conditions may have made long-distance social networks, especially towards Central Europe, particularly prone to distortions. In other words, the phenological patterns probably created large-scale and long-term gradients with effects similar to habitat fragmentation by physical barriers. Importantly, we do not argue for contact vs. non-contact scenarios, where in case of non-contact regional populations become entirely cut-off from larger communication network. Instead, we argue for gradual in- and decreases in contact intensity, the latter fostering increasing variability in the technological and morphological characteristics of the material culture record.

A hypothesis arising from these observations that remains to be tested in the future is that we see an increase in regional idiosyncrasies in the material culture during stadial periods due to more regionalized phenological incentives for populations to move and a decrease during interstadial periods due to intensified contact as a result of more far-reaching phenological pull-factors. However, it must also be taken into account that populations under stadial conditions may have had to move farther to satisfy their needs for subsistence, which would counteract or at least balance weaker phenological incentives for movements.

Conclusion

The presented estimates provide some explanatory potential for spatially and temporally large-scale phenomena regarding probabilities for the distribution and connectedness of Paleolithic populations. Regarding our case studies, we see five main findings.

  1. 1.

    For the expansion of anatomically modern humans into Europe, we see strong phenological incentives for movements along the Mediterranean and Atlantic coasts during interstadial conditions, whereas the route along the Danube shows only weak or no phenological pull-factors. Also, incentives for movements farther into the Iberian and Italian Peninsula are very weak during the interstadials.

  2. 2.

    During the stadial phases, incentives for movements to the southern parts of the Italian and Iberian Peninsula are stronger. At the same time, chains of phenological gradients along the northern Mediterranean coast disintegrate and phenological sinks are very pronounced.

  3. 3.

    Several regions with a particularly dense Aurignacian record, such as Franco-Cantabria, southwestern Germany, and Lower Austria/Moravia/southern Poland, are located in areas where special phenological settings foster a kind of self-contained sub-system, more or less decoupled from the remaining mesh of gradients.

  4. 4.

    In the EEP, radial pulls towards the fringes and north-south-oriented gradients dominate over incentives for east-west movements. In the long run, this makes networks forming in either the western or the eastern part more likely than networks spanning the entire area.

  5. 5.

    Shifts in phenological patterns between interstadial and stadial periods might have been a pace-maker for the growth and disintegration of social networks. In both case studies, interstadial phases seem to be associated with far-flung webs of gradients that foster network growth, while stadial phases are connected with smaller and more compartmentalized patterns that foster network disintegration.

We think that such estimates do not produce arguments for environmentally deterministic scenarios, but rather can provide complementary data and interpretative clues for large-scale temporal and spatial patterns of the Paleolithic record. In the future, we plan to explore this idea further and evaluate the validity of the underlying assumptions by quantifying the net differences in movement incentives and using simulations, such as agent-based models with dynamically changing environmental conditions.