In this paper we present the hitherto first model of potential distribution of Neanderthals across Europe during the Last Interglacial (Eemian, MIS 5e, \(\sim\)130 ka BP). We have applied a robust state-of-the-art modeling approach to take advantage of a sparse set of archaeological records and a palaeoclimatic simulation to produce our model. We analyzed this model using regression-tree based statistical methods to assess how environmental variables shaped habitat suitability at continental and local scales, and assessed to what extent our proposed hypothesis regarding how climate and topography shapes Neanderthal distribution at the continental and global scales matched the results of the model.
The papers by \citet{Richter2005}, \citet{Richter2006}, \citet{Wenzel_2007}, and \citet{Gaudzinski2011} are among the most relevant recent contributions discussing the distribution and adaptations of Homo neanderthalensis. These papers provide a valuable review on Eemian Neanderthal sites (many of them used to calibrate our model), but also promoted a controversial discussion (especially \citet{Gaudzinski2011}) on whether or not the dense forests present in Central Europe during MIS 5e were suitable for Neanderthals due to the lack of large herds of mammals. However, discussion seem to be out of date now, since clear evidences suggest that the Eemian forests of Central Europe were not homogeneously dense \citet{Sandom18032014, Pop201572}. According to our results, the North European Plain, whether it was or not covered by dense forests during the Eemian, was slightly more cold and flat than the optimum Neanderthals’ habitat, resulting in lower habitat suitability values. We do not claim abiotic factors to be more relevant than the structure of the forests in shaping Neanderthals distribution; we are simply adding new information that may help to formulate new questions about this topic. To disentangle the relative effect of climate, topography, and the presence of dense forests across the Neanderthals’ northern range edge would require a reliable land-cover vegetation model, as those introduced by \citet{Allen20102604}, \citet{Gaillard2010483}, or \cite{Fyfe2015}, which are, however, focused on Holocene or Lateglacial time periods.
Two of the papers mentioned above, (\citet{Wenzel_2007} and \citet{Richter2006}) discuss the use of beaches to obtain marine foods (see \cite{Romagnoli2014} and \cite{Romagnoli2015} for examples on the exploitation of marine resources by Neanderthals); more recently, \citet{Cohen201270} also consider coastal areas - many of which are now submerged or have been dramatically altered by the many trans- and regressions of the Late Pleistocene - as highly suitable if not optimal spaces for European hominids. A handful of recently discovered Eemian sites found in Spain and Italy \cite{Arsuaga2012629, Fiorenza_2011}, for instance, are located by the sea, and our model shows that habitat suitability was very high across most parts of the Mediterranean, Atlantic, Adriatic, Aegean, and the Black Sea coasts. In fact, the high habitat suitability values we have found in along the Mediterranean coastline challenges the notion of Neanderthals as a central European species. Our model suggests that the presence records in northern Germany could represent the tail of Neanderthal distribution during MIS 5e, while abiotic conditions were close to optimal along the coastlines. According to our findings, Eemian Neanderthals could have been a species mainly thriving in warm coastline habitats, probably able to exploit the resources provided by the sea \cite{Stringer200814319, Colonese201186} as well as probably small game living near coastal habitats \cite{Blasco2014}. This result is in line with earlier suggestions regarding the importance of coastal habitats in human evolution and dispersal \cite{bailey2004, Cohen201270}. In addition, our rigorous modeling of Neanderthal abiotic tolerance ranges suggests that this hominid species shows no preference for cold conditions.
According to the hypothesis proposed here, the northern edge of the Neanderthals’ distribution should be controlled by cold winters due to the low availability of either small or big game, and cold stress. According to our model, the average temperature of the coldest and the warmest months of the year were two of the three most important predictors at the continental scale. Our local-scale analysis suggests that the lower habitat suitability values are driven by cold winter temperatures below -4.4ºC and positive local-regression coefficients, indicating that these areas (Alps, ID=1; Northern areas of the British Islands, ID=14; Western Norway, ID=40; Scandes, ID=43) are not warm enough during the winter to be suitable for Neanderthals. We further hypothesize that hot and dry summers may have shaped the southern limit of Neanderthals distribution due to heat stress and low summer primary productivity resulting from drought, but also that this effect could have been buffered by proximity to the sea along Mediterranean coastlines. The analysis of our model showed that a few southern localities with Mediterranean influence (Central Anatolia, ID=5; Central plains and South-eastern coast of the Iberian Peninsula, IDs 26 and 27) consistently presented low habitat suitability due to an increased importance of temperature of the warmest month, while habitat suitability was generally very high along Mediterranean coastlines (and the southern coast of the Black Sea).
Finally, our hypothesis puts emphasis on the importance of topographic factors at the local scale: diverse enough to promote a high availability of ecological niches, and hence hunting and gathering resources, but not too rough to fulfill the high mobility requirements of Neanderthals. The analysis of the model at the continental scale reveals that slope is the fourth most important predictor. At the local scale, optimum slope coefficients (<0.029) and slope values between 3 and 11 degrees resulted in high habitat suitability values along the Mediterranean coast, specially in areas with a topographic setting similar to the foothills of larger mountain ranges (e.g. Eastern Pyrenees), calcareous and karst areas (e.g. Karain cave, in the Western Taurus calcareous zone, Turkey) and coastal hilly areas (e.g. Bolomor Cave, in the South-eastern coast of the Iberian Peninsula). On the other hand, the lack of slope (indicated by positive coefficients in the local regression) reduced habitat suitability in plain areas (Pannonian Plain, ID=8; Po Valley, ID=32).
In summary, our model supported the main points of our hypothesis: 1) northern edge was limited by cold winter temperatures; 2) southern edge was shaped by high summer temperature and low water availability (but the latter was only important to separate between high and very high habitat suitability values); 3) high topographic diversity combined with moderate slopes favored occupation at the local scale.
The habitat suitability model presented in this study provides guidelines for interpreting the potential effect of climatic and topographic barriers for Neanderthals dispersal. Mountain ranges like the Alps and Pyrenees probably were strong barriers during the Eemian (and during the early Weischselian, according to \citet{Andel}) that could have offered an opportunity to be crossed during the summer \cite{Richter2006}, but a rough topography could still have been an issue for the mobility requirements of Neanderthals. Extensive flat plains dominated by a continental climatic regime (Pannonian Plain, North European Plain) could also have prevented successful Neanderthal migration, especially during summer in the South and winter in the North.
Our model shows at least six suitable geographical regions, more or less isolated from each other: i) the Atlantic and Mediterranean coasts of the Iberian Peninsula; ii) France and Western Germany, limited by the Alps and the North European Plains; iii) the British Islands, although, interestingly, Neanderthal presence has despite concerted efforts not been confirmed there during the Eemian \cite{Penkman_2011}; iv) Italian peninsula, limited by the Po Valley and the Alps; v) the Peloponese and Aegean Islands, limited by the sea (but see \citet{Richter2005}, \citet{Broodbank_2006} and \citet{Ferentinos_2012}, for a discussion about the possibility of seafaring by Neanderthals), but without any reliable Neanderthal material attributed to MIS 5e; and vi) the coasts of Anatolia, Lebanon, and Israel.
In the long term, the compartmentalization of Neanderthal populations suggested by these relatively isolated distribution pockets could have, especially together with low population densities, led to repeated local extinctions \cite{2013, Premo06012009, PremoB} and genetic differentiation, as shown in \citet{Fabre_2009} where the authors studied genetic variation within twelve samples of mitochondrial DNA obtained from Neanderthal samples dated between 100 and 29 ka BP (but according to \citet{Higham_2014} the latter date may be more probable around 40 ka BP). They defined two genetic demes within our study area (see Fig. 2 in \citet{Fabre_2009}: a northern region that comprises the North of the Iberian Peninsula, Eastern France, Germany and the rest of Central Europe; a southern region comprising the East of the Iberian Peninsula, the South of France, the Italian Peninsula, and the Balkans. These demes do not readily fit with our habitat suitability model: there seems to be connectivity between the eastern and western parts of the Iberian Peninsula due to a corridor of suitable habitat along the south of Pyrenees, and between western and eastern France characterized a large extension of well-connected suitable habitat. But both models agree in the important effect of the Alps as the main geographic barrier shaping the separation between ecologically suitable areas and genetic demes. A detailed analysis of habitat suitability maps along with cultural and genetic evidences can help to understand the causes for gradients of cultural and genetic change across large territories. Likewise, \citet{Ruebens2013341} could identify several culturally distinct regions, each characterized by distinct tools shapes and production methods. As an information transmission system culture is, at the population level, arguably sensitive to forces similar to genetics \cite{Lycett2015}. Like the above genetic analyses, however, Ruebens’ work was also focused on the periods after MIS 5e, albeit with important implications for the interpretation of our abiotic ecological factors.
The modeling approach presented in this paper has certain limitations, most of them grounded to the input data or lack thereof. Neanderthal sites attributed to MIS 5e are scarce in Europe, and their distribution is spatially biased due to a differential availability of areas suitable to preserve such remains (e.g. travertines, caves, lake basins or beach deposits) across Europe \cite{Richter2006}. The chronology of the remains, or more precisely, the method used to define their chronology, is another important source of uncertainty. There are very different methods to assign dates to archaeological remains beyond the range of radiocarbon dating, each associated with its own problems, and usually wide confidence intervals - wider at times than the period in focus here (the Eemian interglacial lasted for around 10,000 years). Therefore, and despite our effort to diminish chronology errors as much as possible by filtering dubious dates, there is the obvious risk of including data points belonging to periods different than the Eemian that may affect the outcome of the model. Weighting the presence records according to their uncertainty may offer one analytical solution to this problem, but one that lies beyond the scope of the present paper. Adding new sites through new fieldwork or improved dating protocols present another solution in the long term.
Our model shows the highest habitat suitability values along coastal areas with mild summers, while the Central European sites, which have dominated our view about the Neanderthal life style during the Eemian, showed low habitat suitability and were interpreted to belong to the distribution tail. Therefore, many current interpretations of Neanderthal livelihood during MIS 5e may not accurately represent this species preferred habitat, i.e. the technological adaptations required to survive in coastal habitats. Further research is needed to investigate the cultural and technological differences between Neanderthal populations inhabiting habitats with such contrasting ecological features.
The methodology presented here to analyze the habitat suitability model is novel, and has proven useful to obtain an insight into the abiotic factors driving habitat suitability at the local scale in different regions across the study area. This methodology can help to avoid discrete representations of presence/absence, using instead SDMs as a source of hypothesis testing and a a source of valuable information for interpreting the ecological conditions at Pleistocene archaeological sites.