Nubian Levallois technology associated with southernmost Neanderthals

Neanderthals occurred widely across north Eurasian landscapes, but between ~ 70 and 50 thousand years ago (ka) they expanded southwards into the Levant, which had previously been inhabited by Homo sapiens. Palaeoanthropological research in the first half of the twentieth century demonstrated alternate occupations of the Levant by Neanderthal and Homo sapiens populations, yet key early findings have largely been overlooked in later studies. Here, we present the results of new examinations of both the fossil and archaeological collections from Shukbah Cave, located in the Palestinian West Bank, presenting new quantitative analyses of a hominin lower first molar and associated stone tool assemblage. The hominin tooth shows clear Neanderthal affinities, making it the southernmost known fossil specimen of this population/species. The associated Middle Palaeolithic stone tool assemblage is dominated by Levallois reduction methods, including the presence of Nubian Levallois points and cores. This is the first direct association between Neanderthals and Nubian Levallois technology, demonstrating that this stone tool technology should not be considered an exclusive marker of Homo sapiens.

. The location of Shukbah Cave, and illustration of the excavations in plan and in section. (a) Regional map of sites in South West Asia and (b) close up map of northern Israeli cave sites preserving fossil Neanderthal specimens, illustrating the position of Shukbah to the south of these sites; (c) a plan of Shukbah Cave (modified from Frumkin et al. 30 ) illustrating areas of excavation, preserved deposits of anthropogenic breccias and the location of the illustrated section shown as (d) (redrawn from Garrod 56  www.nature.com/scientificreports/ In southern Arabia, the presence of Nubian Levallois technology has been argued to reflect to expansions of Homo sapiens from Africa 35,36 . Connections between such technologies and Homo sapiens have also been proposed for the southern Levant 37,38 . Here, we re-examine fossil and archaeological finds from Shukbah Cave, exploring demographic and behavioural diversity in the Middle Palaeolithic of South West Asia. We conduct comparative quantitative analyses of the external and internal structure of the hominin mandibular molar (NHMUK PA EM 3869; abbreviated below as EM 3869) (SI4-6) and Middle Palaeolithic stone tool assemblages that were recovered from the intact, anthropogenic breccias of Shukbah Layer D (SI7).

Results
Non-metric traits of the EM 3869 specimen. The Shukbah tooth (EM 3869) is an isolated fully developed lower right permanent molar (Fig. 2), complete apart from a postmortem chip on the occlusal surface of the distolingual cusp, described in detail in SI4, with the key external and internal features outlined below. Starting with the external description, wear is mild, with pinpoint dentine exposure on the entoconid (wear stage 3 39 ), and there is no distal interproximal facet, meaning that any tooth distal to it had not fully erupted. The regular occlusal morphology, the lack of tapering in the crown shape viewed distally, and the widely spaced and near vertical roots, suggest that EM 3869 is either an M 1 or M 2 . The specimen presents several external (at the outer enamel surface [OES] level, and on the roots) and internal features (on the enamel-dentine junction [EDJ]) that are characteristic of Neanderthals including: (i) the lingual convexity and buccodistal expansion of the occlusal OES outline, in contrast with a straighter lingual margin and buccodistal reduction found in modern humans 40,41 ; (ii) the presence of a wide anterior fovea and a prominent crest linking the mesial cusps (mid-trigonid crest).
The mid-trigonid crest, in particular, is rare or absent in comparative groups for both M 1 and M 2 (SI Table 1) 42 ; (iii) the near rectangular shape of the roots, with little tapering, unlike modern humans, along with a bifurcated mesial root, and buccal and lingual marginal ridges, placed mesially and distally on the mesial root, and mesially on the distal root 43 .
The morphology of the enamel-dentine junction (EDJ) is consistent with the observations of the outer enamel surface, including the presence of a well-developed C5 dentine horn, a deep and elongated anterior fovea and a continuous mid-trigonid crest (SI Table 1 and Fig. 2). EM 3869 also shows an internally tilted metaconid dentine horn, as often seen in Neanderthals 44 , whereas it is more external and vertical in modern humans (Fig. 3a). The relatively short roots and their shape suggest that EM 3869 is an M 1 . The roots of the Shukbah tooth are similar to those of Western Asian Neanderthal M 1 , whilst the M 2 roots are generally longer than the M 1 , less well spaced, and more distally inclined (Tabun I 45 , Amud 1 46 , Shanidar 1 and 2 47 , Kebara II 48 ). We ran a 3D geometric morphometric analysis (3D GMA) of the EDJ shape and tested whether EM 3869 is more likely a M 1 or a M 2 and the results confirm the previous observations, the specimen being statistically attributed to a M 1 (see "Methods"). The presence of an 'X' groove pattern (protoconid and entoconid in contact) at the OES level is both notable and unusual in an M 1 , a 'Y' groove pattern (metaconid and hypoconid in contact) being near ubiquitous in all the comparative samples (SI Table 1). However, Martinón-Torres and colleagues 49 found only 82% with a 'Y' groove pattern in their M 1 Neanderthal sample from largely different sites. At the EDJ level, the morphology of the trigonid crest pattern (type 10) is found in moderate frequency in the M 1 of Sima de los Huesos hominins (21.42%) and Neanderthals (31.25%), but less frequently in modern human M 1 (8.33%) 50 . Overall, EM 3869 exhibits a large part of the morphologically typical Neanderthal features.
Crown dimensions of EM 3869. The size of EM 3869 is particularly remarkable, especially in the mesiodistal length dimension, the adjusted Z-scores for length in all the comparative groups, other than early Homo sapiens, being at 1.0 or above (indicating that the Shukbah specimen is statistically at or outside 95.0% of their  However, many of the teeth from the Western Asian sites are worn, which will have reduced length and increased the figure for crown index. As a subset of Neanderthals, Western Asian Neanderthal M 1 have smaller mean length and breadth dimensions, and the adjusted Z-score of EM 3869 compared with this group for the breadth dimension is at 1.0, in addition to that for the length dimension being at 1.8. Since there are large occlusal wear facets and the entoconid dentine horn tip of EM 3869 is exposed, it suggests that a non-negligible portion of cuspal enamel has been removed by wear. In order to assess tissue proportions, we limited our measurements to the lateral portion of the crown. Lateral enamel thickness estimates are close in Neanderthals and modern humans, the former showing lower average values than the latter (SI Table 7). The Shukbah molar estimates are closer to those of Neanderthals, but also compatible with the range of modern humans (SI Table 7).
Enamel-dentine junction (EDJ) shape of EM 3869. The three-dimensional geometric morphometric analysis (3D GMA) of the M1 EDJ conformation distinguishes between Neanderthals, Pleistocene Homo sapiens, and Holocene Homo sapiens (Fig. 3b, SI Fig. 1 and SI Table 8). The Neanderthal EDJ has a higher topography and more centrally placed dentine horn tips than in modern humans. The specimen EM 3869 falls close to the Neanderthal range, also exhibiting high and internally tilted dentine horns (Fig. 3). Results of the cross-validated bgPCA confirm a high-level of correct classification in the comparative groups and EM 3869 is classified as a Neanderthal (SI Tables 9 and 10).
Root proportions and taurodontism in EM 3869. The roots of EM 3869 are relatively short, display four branches, two mesial and two distal, fused along most of their length. Root length in EM 3869 is a little below the mean for comparative groups of Neanderthals, fossil Homo sapiens and recent Homo sapiens (SI Table 11). The stem root portion goes lower on the buccal aspect (near mid-root length) than on the lingual aspect (only a third of root length). EM 3869 exhibits mild taurodontism (vertical expansion of the pulp chamber into the roots). Taurodontism in other Western Asian Neanderthal M 1 is either absent (Shanidar 47 ; Kebara II 48 ) or mild (Teshik-Tash 54 ; Tabun 1 45 ; Amud 1 46 ). Whilst being associated with Neanderthals, taurodontism also occurs in early Homo sapiens teeth from Skhūl 45 and from Qafzeh, Jebel Irhoud, and Aterian sites 55 . Turning to the internal structure of the specimen, the pulp chamber of EM 3869 shows five well-marked horns corresponding to the main cusps. It is large, expanding low in the root until the bifurcation, before separating into four mesiodistally flattened canals, one in each root branch (Fig. 2). These observations support external assess- www.nature.com/scientificreports/ ment for some degree of taurodontism. We estimated the volume proportions between the stem and branch root portions in the Shukbah molar and it is statistically comparable with the values observed in Neanderthals but differs from fossil and recent modern humans (SI Table 12).
The stone tool assemblage. We analysed 707 lithic artefacts from Shukbah D, constituting 57% of the total collection selected and reported by Garrod 56 and 61% of the collection where its present location is currently known (SI7; SI Table 13), all of which are produced on varying cherts. The results of this analysis support the broad conclusion of Callander's analysis 33 that an emphasis on Levallois point production is evident, but our analysis reveals greater variability within the assemblage than previously acknowledged. Here, we highlight this diversity with respect to patterns of shaping flaking surfaces evident among both the core and flake populations. Our evaluation of the Shukbah assemblage revealed the presence of 12 Nubian Levallois points and 16 Nubian Levallois point cores (Fig. 5). Nubian Levallois artefacts are more numerous at Shukbah than at other sites where they are identified in the Levant and at many sites in Africa 34 , although not as frequent as at some sites in southern Arabia 35,57 . A more exhaustive artefact collection strategy could have further increased their frequency within the Shukbah assemblage. An attribute analysis was conducted to examine morphological differences between Nubian Levallois points and point cores and the wider evidence of Levallois point reduction strategies at Shukbah (see "Methods" and SI7).
Multivariate analysis of the Levallois point dataset indicates that Nubian Levallois points from Shukbah D are not distinct in their overall morphology from the wider body of Levallois points at the site ( Fig. 6a; SI Fig. 2). We extend this analysis and place Shukbah D within its wider context by comparing Levallois point technologies from the site to several key late Middle Palaeolithic comparative assemblages, which date to the timeframe of Neanderthal occupations of South West Asia and/or preserve Neanderthal fossils (SI7). Multivariate analysis of the dataset indicates some inter-assemblage variability amongst these sites, with the notably large collection of Levallois points from Shukbah showing no distinct difference from the total pattern of regional variability ( Fig. 6c; SI Fig. 4). This suggests that selective curation has not skewed the dataset substantially. The results indicate that Nubian Levallois points from Shukbah exhibit morphologies comparable to Levallois points evident across a number of sites dating to MIS 4 and 3 and/or preserving Neanderthal fossils. Comparisons with www.nature.com/scientificreports/ www.nature.com/scientificreports/ southern Arabian assemblages that are notably rich in Nubian Levallois points (123b 58 ; TH383 36 ) as well as Middle Palaeolithic assemblages directly associated with other Homo sapiens (Qafzeh 13 ; Skhul 59 ) indicate significant differences between Arabian Nubian Levallois points and other Levallois points in the sample ( Fig. 6e; SI Fig. 7), which are predominantly explicable through differences in size. The Nubian Levallois points from Shukbah appear morphologically similar to the broader sample of Levallois Points studied, with a subset also overlapping with variability in Nubian Levallois rich assemblages (SI Fig. 8).
Multivariate examination of the core dataset reveals substantive differences between Levallois flake and point cores from Shukbah (SI Fig. 3). The variability of Nubian Levallois cores largely overlaps with the variability observed amongst other point cores (Fig. 6b). Levallois point cores and Nubian Levallois points cores from Shukbah D exhibit comparable variation to point cores present within other later Middle Palaeolithic comparative assemblages (SI Fig. 5), further supporting the potential utility of the assemblage for comparative study. Notably, our inspection of Levallois cores from Bisitun, an undated but excavated cave site in Iran which has yielded a Neanderthal fossil 60,61 , also revealed a Levallois core (Fig. 5i) which matches the description of Nubian Levallois technology 35 , hinting at their more widespread appearance across the region. Further examination of the debitage collections from this may further support this discovery, though in many cases Nubian Levallois technology occurs in low frequency when present 34 . Examination of variability across all Levallois cores sampled from late Middle Palaeolithic comparative sites demonstrates that point cores typically exhibit a subset of variability observed among flake cores, supporting the trend identified at Shukbah (Fig. 6d; SI Fig. 6).
We examined variability among Levallois cores across a range of Late Pleistocene comparative assemblages spanning eastern Africa and South West Asia. Multivariate analysis of this dataset further supports the suggestion that Levallois flake cores are morphologically more diverse than point cores, which show a subset of this diversity, but Nubian Levallois point cores derived from southern Arabian sites demonstrate a departure from this pattern that are predominately driven by differences in core size ( Fig. 6f; SI Fig. 9). Within this context, Nubian Levallois point cores from Shukbah share a common set of characteristics with other Levallois point cores from across the comparative dataset, and extending the variability observed for Nubian Levallois point cores (SI Fig. 10).

Discussion
Our new analyses of fossil and stone tool collections from Shukbah Cave have illuminated unanticipated diversity in both Neanderthal morphology and material culture. Our study provides a thorough evaluation of the EM 3869 molar that enables a robust attribution to a Neanderthal population, corroborating earlier assessment. Based on the degree of wear, the fully developed roots and the lack of a distal interproximal facet, age at death would have been ~ 7-12 years (see SI4-SI5), with wide differences in Neanderthal dental development precluding a more refined age estimate. The size of the Shukbah molar largely exceeds that of other Western Asian Neanderthal specimens from Amud, Dederiyeh, Kebara, Shanidar, Tabun and Teshik-Tash, but more closely approximates the larger crowned European Neanderthals, and in particular those from Krapina (SI Table 3). Indeed, EM 3869 represents one of the largest M1 recovered for this taxon, adding to the morphometric variation amongst Neanderthals. Our firm attribution of EM 3869 as a Neanderthal marks the southernmost fossil of this population known to date.
The stone tool technology from Shukbah D shares the broad characteristics of other later Middle Palaeolithic assemblages found across the Levant, frequently in association with Neanderthal fossils 27,62 . However, our detailed analysis of the Shukbah lithic assemblage has revealed the presence of the Nubian Levallois method, marking the first time this technology has been found associated with Neanderthals. Considerable controversy exists regarding Nubian Levallois technology and its prominent attribution to Homo sapiens 34,63 . The wide distribution of Nubian Levallois technologies across time, space, and technological contexts, as well as typically appearing in low numbers, questions their utility as a robust indicator for cultural inheritance 34,63 . The appearance of Nubian Levallois technologies at Shukbah associated with Neanderthals, and potentially at Bisitun, is most simply explained by technological convergence resulting from a focus on Levallois point production. This is consistent with the broad overlaps seen between Nubian Levallois points and cores and other Levallois point production approaches. Although cultural inheritance of the Nubian Levallois method among Middle Palaeolithic/Middle Stone Age populations cannot be precluded by this finding, they remain to be demonstrated and do not represent the most parsimonious explanation. Our results indicate that any direct link between Nubian Levallois technology and Homo sapiens can no longer be assumed.
The Neanderthal morphology of EM 3869 and the later Middle Palaeolithic character of the stone tool assemblage from Shukbah D are consistent with wider Neanderthal occupations of the Levant between ca. 70-50 ka. Indeed, although Nubian Levallois technology has prominently been associated with MIS 5 by proponents of the 'Nubian Complex' , a detailed evaluation indicates their recurrent appearance amongst Middle Palaeolithic and Middle Stone Age assemblages across the later Middle Pleistocene and Late Pleistocene 34 . As a result, the appearance of Nubian Levallois technology at Shukbah alone offers no clear chronological constraint, and rather the broader character of the stone tool assemblage remains comparable to other assemblages associated with Neanderthals. This is also consistent with the evaluation of the faunal record by Bate 56 , indicating the presence of voles (Cricetinae) and the absence of Afro-Arabian fauna, that is indicative of the wider faunal changes identified at the onset of glacial conditions in the Levant 28,29 . Combined, these support attribution of the occupation of Shukbah D to the wider timeframe of Neanderthal inhabitation of the Levant ca. 70-50 ka, though confirmation of this demands direct chronometric dating.
Examining the limits of the Neanderthal range is critical to evaluate the environmental tolerances of these populations, as well as to explore their ecological adaptability and behavioural flexibility. Our confirmation of Neanderthal occupation at Shukbah presents a notable southward expansion of the range of this population that is associated with diverse technological practices, evident in the presence of the Nubian Levallois reduction www.nature.com/scientificreports/ method amongst other Levallois approaches. These findings help support longstanding assumption of Neanderthal occupations at sites even further south, such as at Tor Faraj and Tor Sabiha 64 . The Neanderthal occupation of the Levant is traditionally associated with Mediterranean woodland habitats, where basecamps appear situated within rough terrain and potentially compress resource niches 22 , but occupation in seemingly arid regions of Syria 65 and Iran 66 suggest wider behavioural flexibility in adapting to South West Asian ecologies. Shukbah presently marks the southernmost point of the Neanderthal range, highlighting the potential importance of this region of the Levant for examining interactions with modern human populations. However, the site may have provided a vital stepping-stone for Neanderthals to expand further south into the arid landscapes beyond. Anthropology System (ASUDAS), and their scoring system has been followed using their reference plaques [67][68][69] .

Methods
The degree of taurodontism of the molar is described using the volumetric bifurcation index 70 Tables 2 and 11). Root stem length was measured as the distance between the cervix and the root furcation buccally. Cervical dimensions were measured as the maximum dimensions at right angles to the mesial and buccal surfaces. The mesiodistal crown dimension of the Shukbah tooth was adjusted for wear using the method of Wood and Abbott 72 . The adjusted length was used to calculate crown area and crown index, and in making comparisons with other teeth. An adjusted Z-score method, using Student's t inverse distribution 73 , was employed to compare each of the Shukbah measurements (external, tissue proportions and volumetric bifurcation index) with the means and standard deviations of comparative groups. The formula applied was: where X, SD and n represent the mean, sample standard deviation and sample size, respectively, of the comparative sample. The interval between − 1 and + 1 comprises 95% of the variation in the comparative sample.
The level of occlusal wear was quantified using Murphy's method, as summarized by Smith 39 . All measurements and observations on the Shukbah tooth were repeated by the same observer after an interval of a month.
X-ray microtomography. The specimen EM 3869 was scanned using the X-ray microfocus instrument (X-µCT) Nikon Metrology HMX ST 225 set at the Natural History Museum of London. Acquisitions were performed according to the following parameters: 183 kV, 195 µA, 3142 images taken over 360° (0.12° of angular step) and 0.7 s exposure time for each projection. The final volumes were reconstructed with a voxel size of 12.75 µm. The microtomographic acquisitions of the comparative fossil and extant hominid specimens were performed using various equipment including X-µCT and synchrotron radiation (SRX-µCT) and reconstructed with voxel sizes ranging from 10 to 42 µm. Data processing. A semi-automatic threshold-based segmentation was carried out in Avizo 8.0 (FEI Visualization Sciences Group) following the half-maximum height method (HMH 74 ) and the region of interest thresholding protocol (ROI-Tb 75 ), taking repeated measurements on different slices of the virtual stack 76 . A volumetric reconstruction was then generated (Fig. 2).

3D lateral crown tissue proportions.
Because of the moderate wear affecting EM 3869, the 3D assessment of tissue proportions was limited to the lateral (non-occlusal) crown portion. A plane parallel to the cervical best-fit plane, and located between the last plane showing a complete ring of enamel and the lowest points of enamel, was used to cut the crown at the cervix. A parallel plane was set at the last point of enamel in the occlusal basin and all material above it (including the part of the crown with the cusps) was removed. The enamel and dentine portions between these two planes were preserved to estimate the lateral crown tissue proportions [77][78][79][80][81] .
The following parameters were then calculated: the lateral average enamel thickness (3D LAET, in mm) and the scale-free 3D lateral relative enamel thickness (3D LRET) 77-81 (SI Table 7). Adjusted Z-score analyses 73,82 were performed on the two tooth crown tissue proportions parameters for EM 3869 in comparison with the fossil and recent hominin specimens/samples presented by Martín-Francés and colleagues 81 (SI Table 7). www.nature.com/scientificreports/  [83][84][85][86] . Five semilandmark curves were set along the marginal outline of the EDJ occlusal basin between each pair or cusps, for a total of 116 semilandmarks. To assess the metameric position, we ran a principal component analysis (PCA), followed by a between-group principal component analyses (bgPCA) based on the Procrustes residuals using two groups (M 1 and M 2 ), as well as a canonical variate analysis (CVA) based on the 13 first principal components explaining ~ 90% of the total variance. Both results of the bgPCA (correct classification: M 1 = 89%, M 2 = 72%; posterior probability for EM 3869: M 1 > M 2 ) and CVA (correct classification: M 1 = 91%, M 2 = 83%; posterior probability for EM 3869: 93% M 1 ) predict EM 3869 is a M 1 . To assess the taxonomic affinities of the Shukbah specimen, we then conducted 3D GMA restricted to the M 1 sample. We computed PCA, followed by a bgPCA based on the Procrustes residuals with the following three groups: Neanderthals, Pleistocene modern humans, and Holocene humans ( Fig. 3b and SI Fig. 1). The bgPCA plot group distribution was verified by computing a cross-validated bgPCA (CV bgPCA; SI Fig. 1 and SI Table 9) and a CVA based on the 11 first principal components explaining ~ 90% of the total variance (SI Table 10). The Shukbah specimen was then projected a posteriori in the bgPCA plot and the posterior probabilities that it belongs to any of the comparative groups were computed with bgPCA and CVA (SI Tables 9 and 10). Results of the bgPCA cross-validation exhibit high classification accuracy, with most misclassifications happening between Pleistocene and Holocene humans, while only one Neanderthal specimen is unclassified (SI Table 9). The CVA give consistent results, with a clear separation between Neanderthals, Pleistocene modern humans, and Holocene humans and the (SI Table 10). Altogether, the analyses are consistent with each other and confirm no spurious and inflated group separation occurs 87,88 . The analyses were performed using the package Morpho v.2.8 89 for R v.4.0.2 90 . Allometry was tested using multiple regressions 91 in which the explanatory variable is the centroid size and the dependent variables are the PC and bgPC scores. No allometry was detected in most analyses (p-value > 0.05), except along bgPC2, where a weak allometric signal (p-value < 0.03; R 2 < 0.1) is detected. The differences between specimens thus mostly represent shape-variation.

Volumetric bifurcation index (VBI) of the roots.
We assessed the degree of root taurodontism in EM 3869 using the 3D bifurcation index 70 . The tooth was virtually bisected into its anatomical crown and root parts by using a best-fit plane at the cervix. An additional plane parallel to this cervical plane was placed at the level of the interradicular surface, dividing the root into the volume of the stem above the bifurcation (Vcervix) and the volume of the branches below the bifurcation (Vbranch). The volumetric bifurcation index (VBI, in %) is calculated as Vcervix/(Vcervix + Vbranch) × 100. We have applied this method to EM 3869 (Vcervix = 502.6 mm 3 ; Vbranch = 234.2 mm 3 ), resulting in a VBI value of 68.2% and compared the results with those of fossil and recent hominin specimens/samples (SI Table 12). Adjusted Z-score analyses 73,82 were performed on the VBI value of EM 3869 in comparison with the fossil and recent hominin specimens/samples presented by Kupczik and colleagues 55 (SI Table 12).  www.nature.com/scientificreports/