Exploring connectivity in Late Bronze Age/Early Iron Age Greece and the Balkans using cranial non-metric analysis

The present study aims to explore connectivity and networking in Late Bronze Age (LBA)/Early Iron Age (EIA) Greece and the Balkans using morphological biodistance analysis and test the potentiality of newly introduced statistical tests, which were designed for challenging datasets, in this particular cultural area. Cranial non-metric traits were recorded in ten skeletal collections, spanning from East Crete to Romania. We followed an experimental statistical approach encompassing two different measures of divergence, the conventional and well-tested mean measure of divergence (MMD) and the newly introduced untransformed measure of divergence (UMD). Though different, results based on these two measures are mutually supporting and show that biodistances in our regional case studies mainly follow the isolation by distance model. This cautiously confirms our main hypothesis that during the LBA and EIA periods in Greece and the Balkans, personal mobility was a slow process characterized by integration, rather than displacement or transformation. The current study is the first one to infer biological affinities using cranial non-metric analysis combined with artifactual evidence, in LBA/EIA Greece and the Balkans. Building a larger dataset through future non-metric analyses will better enable exploring networking and mobility to further complement ongoing bioarchaeological, genetic, and material culture studies.


Introduction
Mobility and migration are themes of regular interest in European and Mediterranean prehistory.The later Bronze Age (1700-1000 BC) is a common focus, when globalizing connections linked groups far and near in unprecedented ways.This study aims to better understand interactions between communities of the western and Southern Balkan peninsula and the Aegean islands.Building on long-standing material culture markers for mobility, through the proxy of human remains, we explore differences between communities during the Bronze Age-Iron Age transition.Connections with the eastern Mediterranean have a long research history, and more recent work has drawn attention to material markers of mobility linking groups in Europe (Borgna 2009;Borgna and Càssola Guida 2009;Bulatović et al. 2021;Iacono 2019;Iacono et al. 2021;Jung and Mehofer 2013;Molloy 2016;Ruppenstein 2020).This is particularly important as the focus of this paper is on communities connected by land or short stretches of water where personal mobility of various modes may be predicted (Anthony 1997).
We studied a sample of human remains in a transect of culturally diverse groups from the Carpathian Basin to Crete.To document how closely cemetery communities were related using human remains, we looked at markers of biodistance.By "cemetery community" we refer to an ancient community represented by their corresponding cemeteries (as an association with settlement remains may not always be possible).The method employed is non-metric trait analysis, which documents morphological variants of the human skeleton (e.g., Buikstra and Ubelaker 1994;Hefner and Linde 2018;Tyrrell 2000).This method provides a point of contact between biomolecular and macroscopic bioarchaeological research (e.g., Herrera et al. 2014;Hubbard et al. 2015;Ricaut et al. 2010).The study is partly exploratory due to the often poor preservation of the human remains we studied, a factor that also impacts detrimentally on genetic datasets also used to evaluate themes of connectedness and biodistance.The human remains studied most commonly date between 1400 and 1100 BC.
In our study area, societal crises around 1200 BC culminated in the collapse of social systems in the Aegean and Carpathian Basin (Molloy 2023).Changes in terrestial and maritime mobility leading up to and following this collapse have been used to model the causes and impacts of this transition in the twelfth to tenth centuries BC (Cline 2021;Emanuel 2017;Hitchcock and Maeir 2014;Jung 2009;Knapp 2021;Knapp and Manning 2016;Knodell 2021;Middleton 2010Middleton , 2017Middleton , 2020;;Molloy 2023;Wallace 2018).Population mobility of various degrees must have been a commonplace phenomenon in the Late Bronze Age Mediterranean world, yet stability and a long-term sense of space and community are also fundamental concepts for the study of this time period (Leppard et al. 2020).
Cultural expression-particularly Mycenaean identitysupported stability and belonging at a political level without implying ethnic or genetic unity."Mycenaean" is an archeological shorthand for distinct population groups that used and exchanged to varying degrees similar material culture, mortuary practices, language, and architecture.There remained a degree of diversity in the hinterlands of societies centered on palaces, and looking progressively farther north, groups selectively used elements of Mycenaean cultural materials alongside non-Mycenaean style objects and practices.We therefore recognize plurality and choice in degrees of "being Mycenaean" within Greece (Barrett and Boyd 2019;Boyd 2016a, b;Feuer 2011).While this cautionary note may be redundant for prehistorians, it is raised due to the conflation of this culture-history term and populations in recent archaeogenetics publications (Lazaridis et al. 2017(Lazaridis et al. , 2022;;Clemente et al. 2021).If being Mycenaean is not implicitly meaningful, in order to understand how different groups within the territory of Greece and the surrounding areas are interrelated, we must evaluate separate communities-through their cemeteries-as building blocks of societies.By focusing on human remains, our intent is to engage with critiques of a mismatch between archaeogenetic and archeological research in studies of personal mobility and culture change (see Gori and Abar 2023 for a recent overview).

The Aegean and South Balkans
In our study, we focus on Late Bronze Age/Early Iron Age Greece and the Balkans.We include a broad cross-section of case studies representing different forms of social organization.Achaea, represented here by the sites of Kallithea 1 and Voudeni, was connected with Italy and parts of non-Mycenaean areas of Greece (Arena 2015;Giannopoulos 2009;Iacono 2019;Jung and Mehofer 2013;Moschos 2009;Tartaron 2014).Although a fortified central structure occurs at Teichos Dymaion, Achaea lacks palatial buildings, and there is as yet no known palatial center (see, however, Paschalidis 2018, 14-15).The community of Kamini Varis in Attika was set within the heartland of the palatial world, bordered by clusters of palaces in Boeotia, Peloponnese, and probably Athens (Osborne 2020).This community occupied the Aegean coast, linking it to the islands, including Crete.East Crete is represented by Halasmenos, Vasiliki-Kamaraki, Istron Kalo Chorio, and Azoria.The Cretan assemblages include remains dating to the LH IIIC to Protogeometric periods and represent the inhabitants of small upland villages or defensible settlements that also participated in long-distance networks (Coulson and Tsipopoulou 1994;Tsipopoulou et al. 2003;Wallace 2018).Monumental tholoi and chamber tombs housing collective burials were dominant burial architecture types during the Mycenaean period in Central and Southern Greece (e.g., Cavanagh and Mee 1998;Boyd 2002, 49-66).Tholoi as well as chamber tombs accommodated collective inhumation burials interred during different funerary episodes (singular and/or multiple, primary and/or secondary burials) (e.g., Moutafi 2021, 30).Similar burial rituals and tomb architecture were used at the Achaean sites, Kamini Varis, and the East Crete group included in this paper.
Other groups farther afield in Epirus and Macedonia in Greece selectively used Mycenaean and Southern Balkan artifacts of material culture as well as a mix of local styles (Bulatović et al. 2021;Eder 2009;Gori and De Angelis 2017;Krapf 2021Krapf , 2017;;Tartaron 2004).These communities occupied a variety of site types, including tells such as Assiros, Toumba, or Kastanas, and enclosed sites, such as Ephyra (Andreou et al. 1996).The community at Achlada was well networked to the north and south (Michael et al. 2021;Ziota 2019).Single inhumations placed in cist/pit graves being the norm in Northern Greece/Southern Balkans contrast with the collective burials common to Central and Southern Greece (e.g., Tartaron 2004).The same pattern of single inhumations placed in individual grave cuts was also found in Achlada (Michael et al. 2021, 26;Ziota 2019).The picture is of a region connected selectively and where material culture use was often highly instrumentalized.Social ties and movement between polities would have followed a range of recognized pathways.

Central Balkans to Carpathian Basin
Looking north to the Central Balkans and the Carpathian Basin, material culture (and pottery in particular) was more extensively used to express social similarity and difference between communities and regions.In seeking to make order of this, the use of culture-history groups has a long tradition.We do not critique this here but rather draw attention to how objects were instrumentalized to express similarity and difference more visibly than in the Aegean area.The southern Carpathian Basin was home to a complex and dense settlement network (Molloy et al. 2020).
A Middle Bronze Age system centered on tell sites gave way between 1600 and 1500 BC to a new social order centered on large open settlements that thrived until ca.1200 BC (Kienlin 2020;Molloy 2023;Nicodemus and O'Shea 2019).Mortuary practices were initially bi-ritual, including both inhumation and cremation, with a bias toward the latter, but this shifted to exclusively cremation in formally organized cemeteries by ca.1400 BC.Given the rarity of inhumations, we treated this area as a unit with individuals coming from multiple sites (Budžak Livade; Pecica site 14; Gakovo-Vasin Do, Rastina, and Mali Akač cemeteries; and the settlement at Turija) (East Pannonian group; see Fig. 1).There was a return to inhumation as a mortuary rite in various forms (articulated and disarticulated) after 1000 BC, including burials at Gomolava and Klisa considered herein.
Along the Adriatic coast of the Balkan region, enclosed forts built of stone are common.Burial in stone tumuli was typical, and remains from Eraci-Ograde, Matkovići, Veliki Vanik, Zavojane-Ravča, Konjsko Polje, Prosik, Zamaslina and Gomila Obličevac (Novak et al. 2011;Lightfoot et al. 2015) are included in this analysis (Dalmatian group, see Fig. 1).In the Balkan interior, tumuli were also the most common mode of burial, and these cover both inhumations and cremations (Dmitrović 2016).While tumuli dominate and cremation extends across the northern parts of this region from the fourteenth century BC, there are some exceptions; these include burials in caves, particularly along the Dalmatian coast and the Lika region, where the Bezdanjača Cave studied herein is located (Martinoia et al. 2021).

LBA/EIA connectivity between the Aegean and the Balkans
Overall, we have a regionally fragmented Bronze Age in Greece, with different communities sustaining connections of variable intensities with each other, parts of Italy, and/or the Balkans and Carpathian Basin for centuries alongside systematic links with communities of the East Mediterranean (Knodell 2021).Exercising choice at a local level was important and so variability in regional identities existed.It is plausible that different regions of Greece had relatively stable populations in which low-level inter-community mobility existed, but directional migration was not a defining element.Leppard et al.'s (2020) macro-scale study of published strontium isotope values for Greece and surrounding areas suggests that communities regularly experienced inward migration of as much as 5.84%.People who had spent parts of their childhood in locations beyond where they were interred would have been present in many nuclear families.Cemetery populations do not correspond in a straightforward manner to living communities in terms of their makeup, but the data suggest regular small-scale mobility covering long periods of residency.Unfortunately due to overlaps in strontium baseline data, it is not possible to say with certainty if this mobility took place within or between Greece, Italy, the Balkans, or the Carpathian Basin (Cavazzuti et al. 2019;Frank et al. 2021;Giblin 2020).An important new resource supporting mobility studies using strontium isotopes is the recently open-access and open-source strontium database SrIsoMed (Nikita et al. 2022), covering countries that have Mediterranean coastlines.
Long-recognized material culture patterns from Southeastern Europe indicate mobility across cultural boundaries, even if we lack precision on the distances, temporality, or scale of such movements (Leppard et al. 2020).Despite a pendulum swing in archaeology back to migration-led narratives, we would argue that this needs to be carefully balanced with a corresponding focus on stability and continuity within communities.Since we struggle to identify instances of rapid or wholesale change in the settlement, mortuary, and material culture evidence, we make use of the consistency of communities in medium-term trends to explore regional population dynamics.

Bioarchaeology to the rescue?
An ideal approach for identifying such medium-scale trends lies in the toolkit of bioarchaeological research, since human skeletal remains constitute the most direct source of information 165 Page 4 of 23 on past populations moving through space and time (Knapp 2021:46).In recent work, radiogenic isotope and ancient genome analyses have more commonly and extensively been used to explore mobility patterns (e.g., Leppard et al. 2020;Marcus et al. 2020;Mathieson et al. 2018;Nafplioti 2008Nafplioti , 2011;;Panagiotopoulou et al. 2018).Even though both methods are powerful paleomobility tools, certain limitations can undermine the validity and representativity of aDNA and/or isotope results.Ethical approaches in the bioarchaeological and biomolecular analysis of human remains are of utmost importance, and destructive sampling in particular requires careful consideration as it is poor practice to sample extensively without preserving relevant remains for future work.This is especially true if we consider the methodological progress over the past 10 years alone, which shows that biomolecular methods are is evolving rapidly and improving to deal with the difficult preservation conditions usually characterizing the south-east Europe and Aegean skeletal assemblages.These conditions pose an additional limitation, namely, the fact that throughout the eastern Mediterranean Bronze Age, the preservation of human skeletal evidence can be uneven or occasionally very poor, and there are still too few palaeogenetic sample sets published to consider regional-scale dynamics.Furthermore, earlier studies on the Mediterranean used mitochondrial DNA (mtDNA) (e.g., Fernández et al. 2014;Hervella et al. 2015) which is neither as robust nor as informative as whole genome analysis, potentially creating sample bias in terms of representativeness (Knapp 2021:48).Nevertheless, aDNA techniques such as whole genome sequencing (Carpenter et al. 2013) and single strand capture (Kapp et al. 2021) have produced promising results.Another limitation of isotopic and aDNA analyses lies in budget capacity, which may result in small sample sizes and consequently in methodological shortcomings (Hakenbeck 2019: 522), such as compromised statistical validity and therefore less meaningful interpretations.A parallel problem is the large spatial scale of certain aDNA projects which may problematically articulate with local-scale studies (Gori and Abar 2023).In this context the analysis of over 30 individuals from Bezdanjača Cave using the capture technique and exhibiting very good DNA preservation is relevant (Lazaridis et al., 2022).The epistemological gap that often exists between archeology and archaeogenetics is also relevant, which is mostly the result of the lack of proper contextualization of genomic data (Ion 2019: 29-30).Finally, other issues that may arise include access to representative samples due to recovery and/or curation conditions.
Macroscopic biodistance analysis on the other hand (e.g., cranial and dental metric and non-metric traits) is less dependent on budgetary challenges as well as access to (and equal collaboration with) laboratory resources and is unaffected by content preservation (aDNA or/and collagen yield).It is also a fully non-destructive approach that returns immediate results.Non-metric traits represent variants of the normal skeletal anatomy that cannot be measured, they are not pathological, and their presence causes no symptoms (Tyrrell 2000).The occurrence of non-metric traits on the human skeleton is commonly documented in bioarchaeology, as over 400 such traits have been identified and they exhibit substantial heterogeneity (Mays 2010).Modeling of non-metric traits has the potential to offer important information on affinities and biological distances between past populations based on the assumption that phenotypic variability expresses phylogenetic variation (e.g., Relethford 2016;Scott 2008).Twin and family studies, as well as studies on laboratory animals, have demonstrated that the majority of these traits are under strong genetic control (Saunders and Popovich 1978;Cheverud andBuikstra 1981a, b, 1982;Carson 2006).Environmental factors also have some influence on their expression (e.g.Grüneberg 1952;Scott and Turner 1997), as supported by the model of quasi-continuous variation (Grüneberg 1952), which became the prevailing explanation for the nature and heritability of morphological traits (Scott 2008).According to this model, the expression of the underlying variation in trait occurrence (genetically and environmentally mediated) is controlled by a threshold; above this threshold, the trait is present and variable in expression, while below it, the character is absent.Under this model, the inheritance of non-metric traits is controlled by multiple genes that interact in order to produce the final phenotype (Nikita 2017, 182).According to Scott and Turner (1997) though, while environmental factors affect non-metric expression, they do not seem to significantly influence the frequency in which the traits occur in a population.What is of interest to our study is that non-metric traits are therefore a combination of inherited and developed features, meaning that they are both genetic and social in their origin.Recent research has demonstrated a correspondence between biodistance results derived from both morphological and genetic data (e.g., Evteen and Movsesian 2016;Irish et al. 2020;Hubbard et al. 2015;Ricaut et al. 2010).Biodistance studies have long been used effectively to differentiate among populations on a global scale (e.g., Hanihara 2013;Hawkey 2004;Irish 1993Irish , 2013;;Turner 1986Turner , 1990;;Scott et al. 2013), as well as on a regional one to address issues such as human mobility and trading networks (e.g., Hubbard 2012;Guatelli-Steinberg et al. 2001;Nikita et al. 2018;Ullinger et al. 2005).
Nevertheless, the aforementioned methodology has limitations.Researchers might face significant difficulties in collecting such traits in a consistent manner, due to preservation issues mentioned above, as well as the funerary practices of the period.Simply put, fragmented remains and incomplete preservation/collection of elements all impact on which traits can be documented.To date, no large-scale non-metric investigation has been conducted in LBA Greece and the Balkans, and a key reason may be due to poor preservation in many assemblages.Mortuary practices are also relevant because the commingling of burials occurred frequently in the LBA Aegean, whereas cremation was the dominant practice in the LBA Balkan-Carpathian area.The comparatively small size of prehistoric cemeteries, the limited number of individuals (minimum number of individuals (MNI)) preserved, and the very large number of missing values in our study area (due to the commonly occurring incompleteness of the skeletal elements) limit the potential of this approach.These traits are also not documented routinely, and even when they are, there is the possibility that different researchers identify or record features in a slightly different manner (inter-observer error).Furthermore, even though the list of cranial non-metric traits proposed by Berry and Berry (1967) is commonly used, it is not as consistently followed as the ASUDAS list of the 21 dental traits proposed by Irish and Guatelli-Steinberg (2003).Therefore, this lack of consistency in terms of cranial traits results in different researchers looking at different traits.We therefore propose the currently used cranial trait selection (see the "Materials and methods" section) as the one that could be more routinely implemented in the specific wider context as a core dataset to build up a larger comparative regional approach.The specific cranial selection offers a large spectrum of alternative possibilities for recording as it includes traits observed on the entire surface of the human skull (e.g., 165 Page 6 of 23 vault, base, maxilla-mandible).This can be most beneficial in cases of assemblages being partially or poorly reserved.It should also be noted that post-cranial non-metric traits were excluded from the current analysis due to the limitations posed by commingling, which very often render the attribution of post-cranial elements to specific individuals impossible.In addition, even though dental traits were recorded, their statistical processing resulted in significant interpretive issues, most likely due to serious limitations (e.g., very small number of dental observations per assemblage; also see the "Discussion" section).Therefore, dental non-metric results will not be presented in the current research paper.
Due to the partial and fragmented nature of most of the skeletons of this context, a large number of cranial nonmetric traits could not be consistently recorded.The main objective of the present study is to explore connectivity and networking in LBA/EIA Greece and the Balkans using cranial non-metric analysis and test the potentiality of newly introduced statistical tests, which were designed for challenging datasets, in this particular cultural area.While the dataset is smaller than desirable, the cranial non-metric traits included were documented in a macro-regional study, and therefore, we present analyses of collections sampled across a geographic scale not commonly applied within bioarchaeological research in this region.We acknowledge that aDNA may reveal exceptional cases, such as migrants from unrelated communities which we can predict were present, and genetic research will provide higher resolution on intracommunity relations, but we are concerned here with the broader background of communities at a more general level.
Taking the state of the art in material culture research as our starting point along with the concept of phenetic distances mainly following the isolation by distance model (Slatkin 1993), two research hypotheses emerge: (1) If immigration into a given community was a gradual process, not fundamentally transformative in the LBA/EIA Greece and the Balkans, the pace of population change would mirror the rate of cultural change, as processes such as integration, accommodation, and adoption-rather than displacementwould characterize human communities.This in turn would allow us to identify distinctive differences between groups of people from the wider region, following the model of low population turnover through inward migration in our period of interest 1400-1100 BC.It would result in biodistances following the isolation by distance model, according to which genetic, and by extension phenetic, distances increase exponentially as geographical distances increase (i.e., positive correlation); (2) in addition, we hypothesize that aside from geographical affinities, communities participating in common cultural routes, based on material culture evidence, would most likely present smaller phenetic distances between them in comparison to culturally more distant or distinct groups.

Materials and methods
This research is part of the European Research Council funded project "The Fall of 1200 B.C.," which investigates migration and conflict at the end of the Bronze Age (ca.1300-1000 BC) and their relevance for understanding collapse in distinct Bronze Age societies of southeast Europe (http:// www.thefa ll1200.eu/ index.html).For the purpose of this particular study, ten (10) assemblages were selected, mostly based on preservation criteria.Cranial non-metric traits were recorded by the first two authors (D.E.M. and L.F.).An inter-observer error investigation would be useful; however, due to COVID-19 restrictions and accompanying very strict regulations in associated museums, this was not possible in this study.Nevertheless, D.E.M and L.F. practiced together on this task (i.e., the recording of non-metric traits) at the University of Edinburgh at the beginning of 2019.Human skulls from the university's reference collection were used to perform this joint training activity, to ensure consistency in practice.Therefore, we believe that possible inter-observer errors are not significantly affecting the results of this study.
The sites that have been included in the current analysis are presented in Table 1 and Fig. 1.Contextual data along with predicted affinities for the ten populations under study are presented in Table 2. Overall, cranial non-metric traits were recorded in 378 individuals from inhumation burials, from 26 sites/assemblages which form ten population groups.Groups were formed under the concept that assemblages that are in very close proximity to each other, contemporaneous in use, and characterized by very similar material culture can be treated as one sample population (e.g., Dalmatian group, East Pannonian group, EIA Serbia group, and East LMIIIC-PG group).A maximum of thirty cranial non-metric traits were recorded for each individual, though in reality, fewer traits were preserved and/or observable in order to be recorded in most cases.The cranial traits were selected on the basis of their successful use by other researchers performing biodistance analyses (e.g., Prowse and Lovell 1996;Hanihara et al. 2003;Sutter and Mertz 2004;Velasco 2018), and the complete criteria for their scoring and dichotomization are summarized in Table 3.Each trait was scored as present/absent, and its frequency was calculated based on the "individual count" method, that is, ignoring potential bilateral presence and only indicating whether the individual presents the trait (Sutter and Mertz 2004).The side with the highest degree of expression was used (Hauser and De Stefano 1989;Scott and Turner 1997).Biological sex was not included in the current analysis, as it would only lead to an additional "fragmentation" of our dataset; therefore, males and females were pooled in the biodistance analysis.Nevertheless, sexual dimorphism was tested as some non-metric Page 7 of 23 165 traits exhibit different frequencies in each sex (Nikita 2017, 188).All age categories (from juveniles to mature adults) were also pooled together, and age dimorphism was explored, as several traits' frequency and degree of expression, such as hyperostotic ones, may increase with age because they involve progressive ossification (Buikstra 1972).Hyperostotic traits are those that are marked by excess bone growth (e.g., tori and bridges) (Saunders and Rainey 2008).Sex and age dimorphisms in trait expression were explored with chisquare analysis (SPSS, version 21).
In order to assess inter-population distances, the mean measure of divergence (MMD) was calculated using the AnthropMMD R package (Santos 2017).Given the fact that in the current analysis, the untransformed measure of divergence (UMD) managed to produce more meaningful results compared to the MMD (see below), all associated information and results on the MMD will be presented in the Supplementary Information.The main difference between the MMD and the UMD lies in the fact that the former is based on a variance-stabilizing transformation formula, contrary to the latter (see Nikita and Nikitas 2021, 39).Due to the mathematical principles underlying its definition (Sjøvold 1977), the UMD measure is an unbiased estimator of population divergence and does not exhibit application problems at very low or very high frequencies (e.g., frequencies less than 0.1 or greater than 0.9), contrary to the MMD (Nikita and Nikitas 2021).Therefore, the standard editing procedure, usually followed in biodistance studies (briefly described in the Supplementary Information), is not necessary when implementing the UMD.The only trait editing strategy consisted of the exclusion of those variables that exhibited only missing values in one or more groups under study.This procedure can be particularly beneficial in challenging datasets, i.e., those with a large number of missing values and a limited number of observations, where most traits would have to be excluded in order to produce an unbiased MMD.Graff 2011;Papadopoulos 1978Papadopoulos , 1980aPapadopoulos , b, 1982Papadopoulos , 1999Papadopoulos , 2000Papadopoulos , 2005 A recent study by Nikita and Nikitas (2021) also showed that the UMD and MMD provided largely the same outcomes and that the UMD generated archeologically meaningful biodistances, in a specific dataset.Nonetheless, this untransformed measure has not been widely tested on large and diverse archeological assemblages and should be treated with caution.In the current study, both measures, i.e., the MMD and the UMD, are used as means of estimating distances and significance, via an experimental step-by-step statistical approach.
Primarily, the biasbinMDs function was used to examine the bias of the calculated UMD measure (for scripts see Nikita and Nikitas 2021).For the stated commands, the Mul-tiOrd, PoisBinNonNor, popbio, and Matrix packages were primarily installed in R. For the estimation of the UMD and its significance the perbinMDs function in R was applied.
Through this method, the accuracy of the p values was estimated.Two statistical tests were used: the S statistic that follows the chi-squared distribution and the T statistic that follows the standard normal distribution (Nikita and Nikitas 2021;Sjøvold 1977;Souza and Houghton 1977).The p values obtained from the stated two statistics were compared to the p values estimated from the permutation method.The permutation method involves the following steps: (1) calculation of the UMD value and the values S and T of the test statistics using the original dataset, (2) random redistribution of all cases into new groups with sizes equal to the original ones, (3) calculation of UMDi and the values Si and Ti of the test statistics based on the new groups, (4) repetition of steps 2 and 3 at least M = 1000 times, and (5) estimation of the p values for the S and T statistics with and without permutations (see Nikita and Nikitas 2021 for more details).Supraorbital notch At least one notch regardless of degree of occlusion Buikstra and Ubelaker (1994: 87) 165 Page 10 of 23 The correlations between the MMD and UMD distance matrices were calculated using the Mantel test (Mantel 1967).Mantel tests have been demonstrated to be robust (Séré et al. 2017), are easy to interpret, and remain widely used to facilitate between-study comparisons (Irish et al. 2020).The same test was used to examine correlations between biodistance and spatial matrices.The latter were based on linear distances.For this measurement, the Mantel function in the Vegan R package was followed, following the script by Nikita (2020, 148).
Phenetic affinities were graphically represented via hierarchical cluster analysis (HCA) and multidimensional scaling (MDS).Regarding the former, dendrograms were used to visualize the relationship between groups and the presence of clusters in the dataset.We also measured the uncertainty in the clustering process, which enables the estimation of the probability of formation of the various clusters and, therefore, their stability (Nikita and Nikitas 2022, 2).The assessment of uncertainty in HCA requires the creation of simulated distances.In order to estimate uncertainty in the dendrograms of HCA, the following steps are followed: first, all pairwise distances among the original samples are estimated, and then, a large dataset of simulated distances is generated; subsequently, the dendrograms of both the original and simulated distances are measured, and the number that each pattern in the original dendrogram appears in the simulated distances is counted.The uncertainty in a pattern is estimated from the percentage of the appearance of this pattern in the simulated dataset, which in turn gives the simulated probability of the formation of this pattern (Nikita and Nikitas 2022, 2).To generate simulated distances, the distance distribution method (DD) was followed as described by Nikita and Nikitas 2021.The Monte Carlo (MC) and bootstrap (B) simulation methods were also implemented to compare with the DD results.All R-scripts and functions are provided in Nikita and Nikitas 2021, while the PearsonDS, mvtnorm, Matrix, and e1071 packages were primarily installed.

Dimorphisms and trait editing
Significant differences were found, for certain variables, between the sexes and between the age categories in trait frequency (Tables 4 and 5).Nonetheless, none of the traits differed consistently in frequency between males and females and between the different age categories in all groups.On the contrary, most of them exhibit significant differences only in one or two out of the ten samples under examination (with the exception of asterionic bone which exhibits significant age differences in three groups).As our dataset already presents certain challenges (limited number of individuals along with a high number of missing values), we decided to primarily conduct the statistical analysis without excluding those traits which show statistically significant sex and age differences, to prevent inflicting further limitations on the biodistance investigation.Nonetheless, MMD analyses were also performed with the exclusion of those traits for comparative reasons.Associated results will be presented in the Supplementary Information.
Regarding pairwise comparison between traits, the statistical analysis showed that only coronal ossicles exhibited a statistically significant inter-correlation.MMD analyses were performed with and without the inclusion of the specific trait, and no differences were observed (MMD values and dendrograms remained the same).

Exploring distances between better-preserved assemblages
As it is obvious from Table 1, most collections of the current study consist of a limited number of individuals, while at the same time, a significant number of missing values are observed in most groups.It should be clarified that the number of individuals stated in Table 1 is the number of individuals where cranial non-metric traits could be recorded based on preservation state and not the overall number of individuals of each site.Bezdanjača Cave, Gomolava, Achlada, Kallithea, Voudeni, and Kamini are the six groups with the highest number of observations (Table 1).

UMD values
When implementing the alternative approach for the estimation of the UMD, which does not require an extensive editing process, 29 traits were selected, with the exclusion of the condylar canal, as it was the only variable that exhibited only missing values in one or more groups under study.The UMD values and significant distances are presented in Table 6.Plots depicting the UMD values are presented in Figs. 2 and 3.As already mentioned, the bias-binMDs function in R was used to assess biases in the estimation of UMD divergence.To test whether a given measure of divergence acts as an unbiased population estimator, the populations from which the dataset originates are simulated.Therefore, the original dataset is used to simulate the respective populations and calculate all pairwise population distances using the UMD (Nikita and Nikitas 2021).Populations of 50,000 individuals in size were created using the φ values, i.e., proportions of the individuals in the samples possessing the trait under examination, as marginal probabilities of the original samples.Then, 1000 samples were randomly drawn from the populations, and their distances were calculated and averaged.The hypothesis was that if a measure of divergence is an unbiased estimator, the averaged sample distances would coincide with the distances between the original samples and between the simulated Page 11 of 23 165 populations; otherwise, the averaged sample distances would deviate from the corresponding original and population distances (Nikita and Nikitas 2021).For this specific case, the measured UMD is a rather unbiased measure of divergence as observed in Fig. 4.
Both MMD matrices (Tables S2 and S3) (with minimum 10 and 22 observations respectively) present positive and significant relations with the UMD matrix (Table 6) (1st case of MMD with UMD: Pearson's correlation = 0.8319, p value = 0.0013, permutations = 9999; 2nd case of MMD with UMD: Pearson's correlation = 0.7544, p value = 0.0013, permutations = 9999).Regarding the correlations between the biodistance matrices and the spatial one, positive correlations were found in all three cases (two MMD matrices and one UMD matrix).However, only the UMD presents a statistically significant correlation with the spatial matrix (Pearson's correlation = 0.4812, p value = 0.044, permutations = 9999).Between the two MMD matrices, the one that was calculated based on more traits and less observations per group (case 1-Table S2) has a higher correlation coefficient with the spatial matrix compared to the MMD of the 2nd case-Table S3 (Pearson's correlation = 0.447, p value = 0.062, permutations = 9999 vs. Pearson's correlation = 0.264, p value = 0.186, permutations = 9999).In terms of networking, both the MMD and the UMD values produced "meaningful" results.In all cases, Kallithea and Voudeni were grouped together, which was expected, based on their proximity in Achaea, Peloponnese.At the same time, the two Achaean groups, which face the Ionian Sea (western axis), present a clear distance and relation with the Kamini group (eastern axis-faces into the Aegean Sea).Gomolava and the Bezdanjača Cave are grouped together in all cases, which again were expected.Achlada in Northern Greece seems to be the "non-stable" element of these comparisons, as it is grouped both with the "Balkan" and Southern "Greek" collections, based on the MMD and UMD values respectively.Interestingly, based on the UMD distances Achlada seems to be "closer" to the Western Achaean sites, compared to the Kamini Aegean site, pointing cautiously to a western axis networking.

Exploring distances between all ten populations
UMD values UMD values were based on 12 traits (metopic suture, supraorbital foramen, supraorbital notch, auditory torus, marginal tubercle, parietal foramina, Inca bone, coronal ossicles, sagittal ossicles, lambdoid ossicles, apical bone, and tympanic dehiscence), as variables that exhibited only missing values in one or more groups under study were excluded.Additional traits with extremely low frequencies were also excluded to produce an unbiased UMD estimator.The UMD values and their significance are presented in Table 7.The graphical representations of the UMD values are depicted in Figs. 5, 6, and 7. Again, it is obvious that the UMD is an unbiased estimator, as depicted in Fig. 8.
The MMD and UMD matrices were positively correlated and statistically significant both when using the Spearman rank coefficient and the Pearson test (Pearson's correlation: 0.7365, p value : 3e − 0.4, permutations = 9999; Spearman's correlation: 0.6102, p value : 0.001, permutations = 9999).Each of the two biodistance matrices presents positive correlations with the spatial matrix, and the results were statistically significant in both cases (MMD and spatial: Pearson's correlation: 0.4474, p value : 0.0087, permutations = 9999; UMD and spatial: Pearson's correlation: 0. 6768, p value = 3e − 0.4, permutations = 9999).In terms of meaningful results, the UMD dendrogram as depicted in Fig. 5 better fits the expected geographical and archeological pattern, compared to the one presented with the uncertainty probabilities (Fig. 7).Both dendrograms employ the same assemblages and depict UMD values; therefore, theoretically, they should have been identical.Nevertheless, they are produced by different R commands and scripts, leading to different clustering.This could be considered a sign of "unstable" groupings.
The four Southern Greek collections were grouped separately from the Northern Greek collections, with Kallithea and Voudeni being grouped together and Kamini being broadly clustered with East Crete (Fig. 5).The latter grouping is quite interesting as both those sites face the Aegean Sea.Furthermore, EIA Serbia is logically grouped with the EIA collection of Gomolava.The two latter collections are broadly clustered with LBA Achlada and the MBA Dalmatian group (Fig. 6) (see Fig. 1).LBA East Pannonian and Bezdanjača Cave groups were broadly grouped with the remaining northern sites (Fig. 5), even though they seem to form distinct groups (Fig. 6).Compared to the MMD groupings, there are some broad and even more specific similarities (e.g., the grouping of Kamini with East Crete); however, it is clear that the extensive although necessary elimination of traits has resulted in the significant loss of clusters and thus information (Fig. S6).
In conclusion, choosing to experiment with the alternative UMD method which does not require an extensive editing process produces more meaningful results compared to the classical MMD approach, especially when including the more challenging collections.It should be noted that the alternative non-editing strategy was also attempted for the MMD; however, the presentation of the relevant results is excluded from this paper since a very large number of traits had to be eliminated to produce an unbiased MMD estimator, as was the case with the AnthropMMD package.

Discussion
According to our original hypotheses, we would expect communities presenting closer geographical proximities and/or participating in common cultural practices to demonstrate smaller phenetic distances between them in comparison to geographically or/and culturally more distant or distinct groups.
The isolation by distance model (Slatkin 1993) was tested through the Mantel test (Mantel 1967).According to Konigsberg (1990), if there is isolation by geographic distance, the correlation between genetic and spatial distances should be positive.In all cases (see the "Results" section), the abovementioned correlations between the MMD/UMD and spatial matrices were positive.Moreover, in most cases, results were also statistically significant.Our results overall do confirm an isolation by distance model, confirming our main argumenti.e., that processes such as integration, accommodation, and adoption-rather than displacement-would characterize community interaction during the LBA/EIA transition in the Aegean and Balkan areas.This allows us to identify distinctive differences between groups of people from the wider region, under the concept of low population turnover in relation to migration.It also shows that blanket terms-such as Mycenaean-used in some genetic studies are problematic, as people using similar material culture were differentiated not only by a number of inherited traits, but by a combination of those traits with a set of specific local cultural routines.For example, Kamini Varis and Achaea, which would both be characterized as Mycenaean, displayed differences in our dataset, reflecting subtle differentiations in their phenotypic variation, also corresponding to variation in their funerary practices (e.g., larger collective burials in Achaea vs. smaller ones in Kamini).This impacts on how we effectively use terminology in an unbiased manner (Eisenmann et al. 2018).Simply put, if our original hypotheses were wrong, the model itself would not work.If there had been high population turnover-for example, even exchanges between those Mycenaean groups in Achaea and Attika-then patterns of biodistance would not correlate with social and/or Euclidian distance.The fact that in few instances, biodistance-spatial correlations were positive but not statistically significant could be interpreted in multiple ways; spatial distances between groups may not be equivalent, with regard to natural barriers (e.g., mountainous areas versus plains); furthermore, gene flow occurred between the communities under study and other populations that were not included in the current analysis.Moreover, the fact that not all populations studied are precisely contemporary could also affect the biodistancespatial results, although exploring connectivity and networking during the LBA/EIA transition lies at the core of our research, and we have been cognisant of material markers for culture change in defining our datasets.Ultimately, this exploratory dataset demonstrates meaningful patterns, but a larger dataset in future non-metric analyses, building on or independent of this one, will better enable the analyses of networking and mobility to complement ongoing genetic, bioarchaeological, and material culture studies.Focusing on our methodology, following a step-by-step "experimental" statistical approach in order to estimate the biodistance measures of divergence (MMD and UMD), MDS and cluster diagrams were produced to visualize phenetic affinities and distances.The potential problem with the latter kind of analysis is that it "forces" samples into clusters; therefore, looking at the MDS "neighborhood" approach might seem more informative.This becomes even more obvious if we consider the quite low simulated probabilities depicted on some of the dendrograms.Nevertheless, we took all information into account when discussing biodistance results.When only the better-preserved six populations are included, the MMD results, in both the examined cases (varying number of traits and observations), agree with the predicted affinities, which are based on geographical location and cultural factors (Table 2).Combining information from both kinds of diagrams, our data modeling suggests that LBA Kamini (Attika) forms a separate group on its own; EIA Gomolava (Serbia) and LBA Bezdanjača Cave (Croatia) form another distinct group or cluster; at the same time, Achlada in Northern Greece lies between the sites of Voudeni and Kallithea in Achaea and between Gomolava and Bezdanjača Cave.When looking at the UMD results, the picture is rather similar with the exception of Achlada that seems to be more clearly grouped with Voudeni and Kallithea.Therefore, it seems that Achlada is the most "fluid" group in terms of biodistances, compared with the other five.Moreover, the UMD uncertainty probabilities for the six better-preserved populations (Fig. 6) are higher compared to those produced for the MMD affinities, possibly suggesting that the UMD is able to produce more stable and potentially more reliable clusters for this specific dataset.Overall, the models are consistent with both cultural contexts and geographic proximity (Table 2).
When including all ten studied collections, the MMD results (when based on only two traits) tenuously produce one broad group including Kamini and the LMIIIC-PG group (East Crete), a broader cluster comprised by the Bezdanjača Cave group (Croatia) and the East Pannonian one (Romanian and Serbian sites) and one tight group with the six remaining assemblages.Therefore, it seems that when datasets are quite challenging, the MMD produces results that are much less informative and reliable.Nonetheless, even in this case, results seem to agree with our predictions (Table 2).In this particular case (ten studied collections), the UMD provided more detailed relationships (Figs. 5 and  6), which are in total agreement with affinities predicted on the basis of material culture analysis.Kamini is grouped with the East Cretan cemeteries, Kallithea is grouped with Voudeni (both sites being in Achaea), and the Bezdanjača Cave site and the East Pannonian group seem to form two separate clusters on their own (however broadly grouped with the other Balkan sites), while four sites, namely, Achlada, EIA Serbia, Gomolava, and the Dalmatian group, form one bigger cluster.The latter is the most interesting and important finding as it could potentially suggest a rather significant phenetic and by extension genetic continuity of communities inhabiting those mountainous Northern Greek/ Balkan sites from the MBA to the EIA periods.Nonetheless, bold assumptions should be avoided, and this result-which will be further tested by our team through genetic analyses in the future-should be treated as merely indicative at this stage.When looking at the UMD and MMD dendrograms with the uncertainty probabilities (Fig. 7 and Fig. S8 respectively), it is obvious that the former seems to be much more consistent with the archeological and geographical context.Nonetheless, its inability to closely cluster EIA Serbia with Gomolava (predicted affinities; Table 2) makes it less meaningful in terms of archeological context compared to the UMD dendrogram presented in Fig. 5. Overall, based on the simulated probabilities depicted in the UMD dendrogram (Fig. 7), it seems that the model is more suitable to distinguish, with much higher stability and reliability, between the locally distanced assemblages, rather than between the more distant ones.
A case that stands out is that of the Achlada cemetery in Northern Greece.As already mentioned, based on different statistical approaches, Achlada lies between the Balkan collections, such as the Dalmatian one, and the two Achaean assemblages in Southern Western Greece.It is tempting to link Achlada's biodistance placement with its geographical location and cultural background.The cemetery of Achlada (near the modern town of Florina) lies south of the historic region of the Pelagonia Valley-an opportune location to trade both north and south, and its location close to the Morava-Vardar-Axios valley is surely relevant (Michael et al. 2021: 22, 26).Combining different kinds of data from its mortuary environment (e.g., artifacts, types of tombs, placement of the dead, adult: non-adult ratio), it can be concluded that features from different traditions were combined in the cemetery (see Michael et al. 2021 for more details).The material culture shares elements with finds from cemeteries in Ulanci, Central Macedonia, and Albania as well as objects reflecting local traditions (Michael et al. 2021;Ziota 2019).Molloy (2016, 369) argues that the material evidence for land-based interactions in the LBA Southern Balkans points to combinations of regular and sustained links with neighboring groups and eclectic or irregular encounters with more distant groups.The former may have included marriage exchanges, for example, whereas the latter may have included both circular (return-to-home) and chain migration, to use Anthony's (1997) terminology.A potentially east-west axis of influence may be evident between the Aegean-Ionian Seas, as suggested by the current analyses.Kamini in Attika is clustered (in a broad sense) with East Crete, and this pairing seems consistent in most dendrograms.The Attika assemblage is located on the Aegean coast, linking this community with the islands, including Crete.Kallithea is placed closer to some assemblages to the north and northwest (such as Achlada or the Dalmatian group) than to Kamini Varis in Attika (e.g., Fig. 6; Tables 6 and 7).A similar pattern to that seen for Kallithea is observed for Voudeni, also in Achaea, but only in the UMD values (see Table 7).Achaea was connected with Italy and non-Mycenaean areas of Greece to the north (e.g., Iacono 2020).Therefore, it seems archeologically reasonable that when compared with the specific groups of our research, the Achaean sites cluster closer to western than eastern datasets.Nonetheless, larger and broader assemblages should be included in future non-metric analyses to form a clearer picture regarding interactions between LBA Greece, the Balkans, and potentially Italy.
Recent genetic work on LBA individuals from Southern Greece and Crete has identified spatially relevant internal variation, and work in Crete has defined temporally relevant changes in population genetics (Clemente et al. 2021;Lazaridis et al. 2017Lazaridis et al. , 2022;;Skourtanioti et al. 2023).Furthermore, Lazaridis et al. (2022) argue that several individuals from Bulgaria had "the Mycenaean genetic profile, which suggests that Mycenaeans were genetically similar to" some contemporary people there.Though Lazaridis et al. unhelpfully speak of "Mycenaean Greeks" as a meaningful entity without defining spatial boundaries, it is clear that variation must increase as one moves north since LBA Carpathian populations are shown to have a defined genetic distance.However, there is a bias in sampling creating a dislocation between groups studied in the south and the north.LBA samples come from the far south of Greece and from the central to north Carpathian Basin, with virtually no samples included from parts of the Balkans lying in between.Despite this bias, Lazaridis et al. (2017Lazaridis et al. ( , 2022) ) emphasize genetic difference by capitalizing on this knowledge gap to support a narrative that takes insufficient account of archeological literature and has been robustly critiqued (Frieman and Hofmann 2019;Gori and Abar 2023;Hakenbeck 2019;Hamilakis 2017).Given the predictable diversity of groups across this tract of Southeastern Europe, this ill-defined disjuncture limits our view of long-distance interactions, mobility, and relationships.More specifically, it restricts analysis of the relative role of these two centers of cultural significance-the Aegean and the Carpathian Basin-on each other and on the lands in between.The current study is a step toward exploring the breaks in this chain and will be followed up by further genetic analyses in preparation by the authors and collaborators at the Globe Institute, Copenhagen.
Our results support a model of integration of migrants rather than rapid transformation or displacement, which aligns well with recent discussions on mobility and migration in archeological studies (Heyd 2017;Iacono 2019;Knappett and Kiriatzi 2017).This picture also fits well with Leppard et al.'s (2020) isotope analysis paper, estimating 5.84% of non-local presence in the ancient Mediterranean between the MBA and EIA (2000-800 BC).Leppard and colleagues, reflecting on the arguments of Iacono ( 2019) and Legarra Herrero (2016), concluded that if the whole idea of a connected Mediterranean was indeed a reality, it was limited to those of the higher social strata, whereas for the average person, long-distance or frequent mobility, and certainly migration, was exceptional (Leppard et al. 2020, 231).Looking at our contextual data (Table 2), of our individuals stand out as striking, which may suggest that we are looking at typical rather than exceptional members of a community.On the other hand, certain limitations in the strontium analysis methodology should not be overlooked (e.g., overlapping strontium values between regional zones; the possibility that the founding generation of migrants is not being well represented; underrepresentation of datasets in specific areas; see Leppard et al. 2020, 227-230).Should we accept the 5.84% mean non-local persons model that Leppard et al. (2020) propose, this translates roughly to one out of every 20 individuals in each community being a migrant.Even if we assume that these communities were primarily composed of small family units of four people (two adults and two children), one individual in every five families may be a migrant, allowing for larger family units, that are more likely to be one in every two or three families, which could support the concept of small-scale mobility on a regular basis.
In sum, cranial non-metric analysis was able to provide insights on networking and mobility in this cultural area.The current study has limitations, including a restricted number of observations and a large number of missing values, mostly due to preservation conditions.Such limitations are particularly evident in specific dendrograms, which present very low simulated probabilities (e.g., Fig. S8).Nonetheless, even in such cases, the model is able to produce archeological meaningful pairs (e.g., Kallithea and Voudeni).Dataset limitations also led to the inability of incorporating the factor of biological sex or even age into the research design.Statistically significant sex and age differences in trait frequency were found; however, none of the traits differed consistently in frequency between males and females and between the different age categories in all groups.Even though excluding such traits from the analysis led to further challenges (see the "Results" section) in certain cases, results were identical to the ones produced without eliminating those traits, possibly demonstrating a certain level of stability.It is evident, however, that apart from securing a minimum of 10 observations per group (Irish 2016), incorporating as many traits as possible into the non-metric analysis produces more trustworthy and archeologically meaningful results.
The current research used both the well-tested MMD as well as the newly introduced UMD to infer biological affinities and distances.Both approaches managed to produce archeologically meaningful results; however, the UMD approach led to more detailed evidence of relationships between the ten groups under study as well as higher simulated probabilities, possibly suggesting more stable clusters.Even though the UMD succeeded in producing meaningful results in both our study and the one of Nikita and Nikitas (2021), it has not been tested independently by other researchers with expertise in statistics, and more importantly, it has not been employed in larger and more diverse assemblages which form known biodistance clusters (Nikita and Nikitas 2021, 40).It should also be noted that when we tried to explore biodistance through dental non-metric traits, it was not possible to produce archeologically meaningful results with either of the two employed strategies (MMD and UMD).This was possibly due to the inclusion of even more problematic datasets (for example limitations of associating mandibles with maxillae in commingled deposits) and the recording of fewer traits compared to the number of cranial ones recorded for this research (21 dental traits vs. 30 cranial).Therefore, even though the current study primarily supports the alternative strategy as outlined in Nikita and Nikitas (2021), it also stresses the importance of investigating each dataset carefully and recognizing its unique characteristics.Our future goals include integrating non-metric biodistance data into isotope (both stable and radiogenic) and paleogenomic analyses, to form a holistic picture of connectivity in LBA/EIA Greece and the Balkans.Regarding the latter type of analyses, current conclusions will be tested through genetic investigation of the same cemetery collections undertaken by Hannes Schroeder and Miren Iraeta at the Globe Institute in Copenhagen.

Conclusions
Biomolecular research will continue to become more affordable and grow as a powerful field for exploring biodistance.Yet its destructive nature has much potential to pit it against ethical and sustainable research values and codes developed over decades of bioarchaeological research.Taking a different but complementary approach to biodistance, this study is the first to infer biological affinities using cranial non-metric analysis, combined and compared against studies of artifactual evidence, in LBA/EIA Greece and the Balkans.Using 10 skeletal assemblages, spanning from East Crete to Romania, this study provides provisional insights on networking and mobility on a local and regional scale while properly acknowledging inherent limitations.We followed an experimental statistical approach encompassing two different measures of divergence, the conventional and well-tested MMD and the newly introduced UMD.Uncertainty probabilities in HCA dendrograms were also estimated.This study supports the UMD approach in biodistance studies and stresses the importance of further exploring this alternative measure of divergence in larger datasets with known affinity clusters.Though different, results based on these two measures of divergence are mutually supporting in the manner argued above.These show that biodistances in our regional case study mainly follow the isolation by distance model, which in turn cautiously confirms our initial hypothesis that during the LBA and EIA periods in Greece and the Balkans, the contribution of mobility to changes within a community was a gradual process characterized by integration, than displacement or transformation.This is in alignment with archeological evidence but also with insights into mobility patterns revealed by strontium isotope studies in the wider area as well as recent genetic work.Building a larger dataset through future non-metric analyses will better enable the exploration of networking and mobility to further complement ongoing bioarchaeological, genetic, and material culture studies.

Fig. 1
Fig. 1 Sites of modern Greece and the Balkans that have been included in this study

Fig. 2
Fig. 2 Ward's dendrogram depicting the UMD values also including the uncertainty probabilities (in blue font); (the same dendrogram was produced following the perbinMDs function as described in Nikita and Nikitas 2021)

Fig. 4
Fig. 4 UMD values (following Nikita and Nikitas 2021) between the averaged sample distances (mean), the distances between the original samples (sample), and those between the simulated populations (pop), for the six better-preserved groups.The averaged sample distances coincide with the distances between the original samples and between the simulated populations

Fig. 7
Fig. 7 Ward's dendrogram for 10 groups based on the distance distribution method following Nikita and Nikitas 2022 and including the uncertainty probabilities (in blue font)

Fig. 8
Fig.8UMD values between the averaged sample distances (mean), the distances between the original samples (sample), and those between the simulated populations (pop) for the 10 groups.The averaged sample distances coincide with the distances between the original samples and between the simulated populations

Table 3
Cranial non-metric traits recorded for this study

Table 4
Significant sex differences in cranial non-metric trait frequencies