On the Japanese archipelago, the Jōmon period (14,000-800 cal BC; Incipient: 14,000-9300 cal BC; Initial: 9300-5200 cal BC; Early: 5200-3400 cal BC; Middle: 3400-2400 cal BC; Late: 2400-1200 cal BC; Final: 1200-800 cal BC [1]; the Incipient period was excluded in this research because no human skeletal remains from this period were found), when people engaged mainly in hunting and gathering, is prehistoric period that is particularly important for researchers investigating how people reacted to environmental fluctuations in human evolutionary history. This period experienced more severe climatic changes, including the Younger Dryas, than did subsequent Yayoi (800‐250 AD) and Kofun (250‐600 AD) societies [2-5]. Many have claimed that prehistoric societies undergoing radical environmental fluctuations tended to be accompanied by societal change [6-9]. The Jōmon society may be no exception. Research has suggested that climate change has affected population size in many regions [10,11], which may have influenced other aspects of the Jōmon society.

However, features such as the emergence of social hierarchies, shifts in subsistence patterns, and significant technological developments commonly seen in other periods were less radical in the Jōmon society. Jōmon people predominantly adhered to a hunter-gatherer lifestyle and persisted in using stone tools for over 10 millennia. Although recent findings at the Usu-Moshiri site in Hokkaido found evidence of violent injuries in more than 5 individuals, such occurrences remained relatively rare within the broader context of the Jōmon period. This stands in stark contrast to the following Yayoi period [5,12-16] and can be distinguished from patterns observed in other relevant regions worldwide [17-19]. Although the notion of a “static view” of the Jōmon society has been challenged over the years [2,20], it remains clear that societal changes during this era were relatively less drastic, even in the face of more drastic climatic variations.

One potential strategy for mitigating the impact of external environmental stresses is sustained interaction between different populations. The continuous exchange of genetic material, knowledge, skills, and resources can enhance the resilience and adaptability of populations in the face of environmental disruptions. Many relevant studies of genetic diversity in organisms other than humans support this possibility, including plants, various animal species, and whole ecosystems [21-23]. The specific mechanisms behind the genetic diversity and resilience remain controversial, although it is generally argued that genetic diversity is somehow related to morphological and functional diversity, and that damage to populations due to ecological disturbance could be compensated for by such diversity.

This is also evident in modern hunter-gatherer societies, which exhibit higher genetic diversity, indicative of ongoing interactions and gene flow between groups [24-26]. Moreover, ethnographic research shows that mobile hunter-gatherer communities with lower population densities often seek mates from distant groups [27-29]. Archeological and anthropological research also suggest that higher mobility of prehistoric hunter-gatherers was an important factor for their societal stability [11,30,31].

Anthropologists have extensively discussed the population interaction history during the Jōmon period of the Japanese archipelago. Typically, previous research used a traditional biodistance method using distances between specific measuring points[32,33] (see also Figure 1), and statistically significant differences between morphological biodistances among different populations have been regarded as populational and genetic differences. For example, as stated in the next paragraph, if geographical clines and statistically significant differences are found among morphological biodistances of different populations, it is commonly suggested that populations moved in a certain geographical direction and that they did not migrate widely.

Some research on human skeletal remains including crania and limbs has supported the notion of relative morphological uniformity within the Jōmon population, which could be sustained by widespread population interactions, although others have claimed that spatiotemporal differences may exist to some extent [34-42]. Dodo [34,35] argued that statistical differences found in cranial biodistances among Jōmon populations in the Tohoku, Kanto, and western regions were much smaller than among Jōmon populations and more recent populations. On the other hand, it has often been claimed that geographical clines are found from north to south because the Jōmon people were assumed to be influenced by both the southern and northern populations outside the Japanese archipelago. Fukase et al [36] claimed that such geographical clines were found in the length of limbs. Hanihara and Ishida [37] and Kondo [38] argued that geographical clines were also found in cranial biodistances. Ogata [40] and Yamaguchi [41,42] showed that the body size of the Incipient and Early Jōmon people was smaller than that of the other subperiods due to the dietary conditions.

Many archeological studies have also discussed population interactions via a widespread distribution and exchange of artifacts such as pottery and lithics over extensive geographic regions during this period [2,3,43-45]. It has often been pointed out that a certain type of pottery was widely distributed across regions. An example is the Ento-kaso type distributed in both the Hokuriku and the northern Tohoku regions of the early phase (see Figure 2 for regions in the Japanese archipelago). Even if the distribution of a type of pottery was regionally restricted, it may have been affected by a distant region. The Funamoto type of the Middle phase was distributed mainly in the Sanyo regions, and this type of pottery strongly influenced the pottery in the Kyushu region during that same phase. In the Late phase, when samples were obtained from the largest number of different regions in this study, it was found that a large quantity of the Horinouchi type I pottery originating from the Kanto region was excavated in the Kinki region; similarly, the Nakatsu type pottery from the Kinki region was transported to the Kanto region and influenced the origin of Shomyoji pottery [45,46]. Such diverse and widespread distribution of pottery suggests that the Jōmon people interacted widely, although it has often been argued that Jōmon society was regionally diverse and distinct, with differences in settlement styles, lifestyles, and types of archeological remains, such as clay figurines. This is why the Jōmon cultural complex has been referred as “Jōmon cultures” [1,47,48].

The geometric morphometrics of 3D data, the approach used in this research, which explores shape variations of targeted objects typically using the coordination of landmarks, has gained traction in various fields, extending to objects such as lithics [49-51]. Although the morphological variations of the Jōmon people have been explored in previous studies using the traditional biodistance method, this study has the following novelties. The 3D data of Jōmon crania have rarely been examined using geometric morphometrics, with some exceptions [52-56], and these exceptions focused on smaller samples from restricted regions. Less biased larger samples are important for estimating macroscopic population interactions. Furthermore, when analyzing the resulting data, geometric morphometrics is able to track morphological variation among configurations of each measuring point at once, whereas the traditional biodistance method compares each measured distance independently. These differences could provide some new insights.

This study aims to examine the hypothesis that the Jōmon populations interacted widely and continuously in an anthropological way (ie, through geometric-morphological analysis of the cranial features of the Jōmon people) using a much larger 3D cranial sample collected from broader regions than previous studies. If the populations interacted widely and continuously, there should be fewer morphological differences among different regions and phases (ie, fewer statistically significant differences between each region and phase are found). Moreover, when the Jōmon populations are compared with populations from a different period, the hypothesis would expect that interregion and interphase variations should be lower than in the populations from a different period.

We chose 31 representative and/or less complicated landmarks for this geometric morphometric study, which were placed on 3D models using the geomorph package in R (version 4.3.2; R Foundation for Statistical Computing; Figure 1). Our research focused on the human cranium, enabling us to use the landmark method, which was used in similar studies [54-58]. Landmark selection was based on those described by Martin [32] and Caple and Stephan [33] (see Figure 1 for the Bookstein type of landmark [33] selected in this study).

We collected relatively well-preserved and less skewed 3D cranial data from the Jōmon period in collaboration with local archeological centers, museums, and universities. We used a range of scanning technologies, specifically Creaform HandySCAN BLACK, HandySCAN BLACK Elite, and Einscan HD, in addition to Structure from Motion and Multi-View Stereo techniques based on 2D photographs [59,60]. Previous studies have confirmed the consistency of the models generated by the aforementioned methods for human crania [60,61]. Our dataset comprised 363 Jōmon crania (including 146 females, 215 males, and 2 individuals with unknown sex), with specific information available in Table 1 and Multimedia Appendix 1. These 3D data were collected from various regions spanning the island of Honshu within the Japanese archipelago, encompassing regions from Tohoku to Kyushu (Figure 2). Although buried crania from the prehistoric period tend to be distorted due to soil conditions, we chose to measure crania that are regarded as less skewed and better preserved based on their appearance and relevant excavation information.

To provide a basis for comparison, we included cranial data from the Kuma-Nishioda site of the Middle Yayoi period (350 cal BC to AD 30) in the Kyushu region and the Aoyakamijichi site of the Late Yayoi period (AD 30 to 250) in the Sanin region, comprising 22 and 18 crania, respectively. The choice of the Kuma-Nishioda site was motivated by the common assumption that Yayoi people, especially in the northern Kyushu region, might exhibit lower genetic and morphological diversity due to a genetic bottleneck caused by gene flow from migrating Korean populations [55,56,62,63]. The Aoyakamijichi site was selected because it is the only one with a relatively large number of 3D samples outside of the Kyushu region [64]. Since we could not gather good 3D data from the Early Yayoi period, we did not include this period.

The file size of the original 3D data was so large that it took a long time to read the data and locate the landmarks in R; the files ranged from 100 to 600 MB depending on the measurement methods. To address this issue, we reduced the size to 2.5%‐15% of the original using MeshLab [65]. Previous studies have confirmed that this reduction process does not result in a model that is significantly different from the original [59,66,67].

Some prehistoric crania were incomplete, so we attempted to reconstruct lost areas in 3D models. This reconstruction involved mirroring remaining portions along the median line, passing through the nasion, prosthion, and bregma [32,55,56,68,69]. Since this mirroring could include some subjective biases, mirroring was limited to cases with smaller distortions and clearer original positions of remaining parts. In cases where mirroring was not possible, particularly when corresponding parts were entirely lost, we applied the estimate missing function, using the thin-plate spline method within the geomorph package (version 4.0.4) of R (version 4.2.1; R Foundation for Statistical Computing) and R studio (2022.07.0 + 548) [70-72].

Sex and age estimation of skeletal remains depended mainly on published excavation reports, where anthropologists estimated the sex and age. If any relevant information was not given in the reports, we made educated guesses by investigating pelves and crania (especially postcranial features and sizes of the zygomatic arch and mastoid process) for sex [73], while the pubic symphyses [74] and cranial sutures were used to estimate age [75,76].

After conducting generalized Procrustes analysis to rotate and adjust the landmark configurations for all Jōmon and Yayoi crania, we performed principal component analysis (PCA) to investigate the spatiotemporal variations of these crania. We did not perform PCA for each phase or region because no region had samples from all phases, and no phase had samples from all regions. Additionally, we also performed the Steel-Dwass test, which is a nonparametric multiple comparison test with an expansion of the Mann-Whitney U test, to evaluate the principal component (PC) scores of Jōmon crania based on each phase and region. All statistical analyses were performed using R and R Studio.

We also conducted the PCA and the same statistical tests on all the Jōmon and Yayoi populations to clarify how varied the Jōmon populations were when compared to populations from a different period. Next, we selected cranial data from Jōmon sites with relatively larger sample sizes to compare with the Kuma-Nishioda site. They were compared via PCA by each site. The latter analysis included the Ebishima (from the Tohoku region of the Middle phase), Nakazuma (from the Kanto region of the Late phase), Yoshigo (from the Tokai region of the Final phase), and Tsukumo sites (from the Sanyo region of the Final phase), and their locations are indicated in Figure 2 and Multimedia Appendix 1. Such comparisons of the Jōmon populations from different phases and regions with the Yayoi population could reveal the degree of variation of individual Jōmon populations.

Ethical approval was not required for this study according to the local legislation and institutional requirements because under the Law for the Protection of Cultural Properties in Japan, human skeletal remains are regarded as cultural properties [77]. Curators of all samples are described in Multimedia Appendix 1.

The PCA results indicated cumulative contribution rates exceeding 75% up to the 20th PC (Table 2). To streamline the discussion, this study primarily focused on PCs with magnitudes larger than 5% (PC1-PC5). PC1 was positively correlated with narrow facial width, prognathism (protrusion of the jaws), and high and short occipital areas, while PC2 was positively linked to temporal length and negatively related to facial length. PC3 primarily captured temporal length, and PC4 reflected facial height. PC5, primarily representing cranial asymmetry, was not included in the analysis as it was influenced by ground pressure during burial (Figure 3). We also showed deformation patterns of a cranium from the Fujizuka shell midden (ID 341; see Multimedia Appendix 1) according to the morphological change each PC captured (Figure 4).

The PCA results, graphically presented in Figures 5 and 6, suggested that the differences between phases are relatively small, in line with the statistical tests on PCs that indicated statistical differences only in the Early and Final (Z=3.118, P=.02), Middle and Final (Z=4.233, P<.001), Late and Final (Z=4.040, P=.001) in PC2, and Late and Final phase (Z=2.946, P=.038) in PC4 (see Multimedia Appendix 2). Notably, the scatter plots in Figure 5 and the box plots in Figure 6 showed more prominent regional differences, consistent with the statistical tests on PCs showing significant distinctions between various regions. These differences were particularly notable between the Tohoku and Tokai regions (Z=6.375, P<.001), the Tohoku and Sanyo regions (Z=5.852 P<.001), the Kanto and Tokai regions (Z=4.880, P<.001), the Kanto and Sanyo regions (Z=4.180, P<.001), the Hokuriku and Tokai regions (Z=3.773, P=.004), and the Hokuriku and Sanyo regions (Z=3.829, P=.003) in PC2. Additionally, there were significant differences in PC4 between the Kanto and Shikoku regions (Z=3.302, P=.02) and the Shikoku and Kyushu regions (Z=3.088, P=.04). It should be noted that there were no significant differences observed in PC1 and PC3, and the regional differences remained relatively limited. Furthermore, there were few sexual differences (see Multimedia Appendix 3), with some exceptions, specifically PC3 and PC4 in the early phase of the Tokai region and PC4 in the middle phase of the Kanto region, which was supported by statistical tests between the sexes in each PC score (U=15,926, P=.81 in PC1, U=13,925, P=.07 in PC2, U=12,693, P=.002 in PC3, and U=13,695, P=.004 in PC4 comparisons).

The two comparative results between the Jōmon and Yayoi populations were summarized in Figures 7 and 8 and Multimedia Appendices 2 and 4. In the second site-to-site comparison, to streamline the presentation, these figures showed only PC1 and PC2, as all sites displayed substantial overlap in PC3 and PC4 except the results between the Yoshigo and Kuma-Nishioda sites (see Multimedia Appendix 4).

The statistical results of the first comparison indicated that statistically significant results were found only in PC2 and PC3. PC2 focused on the interphase and interregion differences among all populations, and PC2 and PC3 captured the difference between the Jōmon and Yayoi periods (Multimedia Appendix 2). Figure 7 showed that the Jōmon populations substantially overlapped, while the Yayoi populations were more scattered spatiotemporally in PC2 and PC3, which was supported by the fact that the effect sizes of statistical tests between the Yayoi populations (Z/n=.113 in spatial and temporal comparisons in PC2) were higher than the ones between the Jōmon populations (Z/n<.055 in spatial comparisons and <.029 in temporal comparisons in PC2; Multimedia Appendix 2).

Figure 8 showed that the most pronounced differences of the second comparison were observed in PC1, with one notable exception seen when comparing the Ebishima and Kuma-Nishioda sites, where the differences were more conspicuous in PC2. Overall, the Jōmon individual populations were more varied than the Yayoi population. Figure 9 shows a visual representation of the configurational changes in landmarks for each comparison. The figure was constructed using the plotRefToTarget function of the geomorph package in R. Straight bars indicate the degree of variation in each landmark. These overall differences are primarily related to facial height and tooth position, specifically highlighting that individuals from the Kuma-Nishioda site tend to exhibit taller facial features and a more anterior tooth placement.

The results of the PCA presented herein were not straightforward. PC2 showed statistically significant spatiotemporal differences between some regions and phases, which was consistent with previous biometric studies claiming interphase variation [40-42] and geographical clines from north to south [36-38,41], although it is difficult to conclude that their differences are pronounced. This is because PC1 did not show any statistical differences, and a similar study of the Kofun period, examining a larger set of 3D data of human crania using geometric morphometrics, exhibited a geological cline in PC1 [56]. Moreover, as mentioned in the Introduction, given the archeological evidence for the spatiotemporal distinctiveness of the Jōmon material cultures, it would be expected to find statistically significant differences in PC1 as well. Thus, it is possible that the spatiotemporal differences throughout the Jōmon period were more nuanced or relatively small.

Comparisons between the Jōmon and Yayoi populations are consistent with the above interpretations. It is suggested that morphological variations among Jōmon populations are relatively less varied spatiotemporally than in the Yayoi populations, and individual populations are more diverse than the Yayoi population from the Kuma-Nishioda site. It has been widely proposed that Yayoi populations, particularly in the Kyushu region, are descended from migratory groups originating from the Korean peninsula [55,62,63]. Suppose that a genetic bottleneck resulting from this migration contributed to the reduced variation observed in the Yayoi people of northern Kyushu. It suggests that the morphological and genetic diversity among the Jōmon populations were not relatively limited.

Examining population interactions in the Jōmon period could provide important insights into their resilience to environmental fluctuations at the time. As stated in the Introduction, the Jōmon period experienced more severe climate changes than subsequent periods such as the Yayoi and Kofun periods, which sometimes showed a radical societal change. The reason why less evidence of societal changes—such as the occurrence of warfare and emergence of social hierarchies—was found in the Jōmon period may be because the Jōmon populations interacted widely and continuously and frequently exchanged knowledge, skills, and resources against environmental disruptions.

First, even if our results suggest that Jōmon people moved or interacted widely, it is still possible that such movements or interactions were not continuous but rather discontinuous with some drastic environmental changes such as the Kikai-Akahoya eruption and the 4.2 ka event [78-80]. The frequency and degree of interactions naturally depended on the region and phase, as suggested by research on pottery distributions [45]. It should be noted, however, that certain regional interactions must be maintained to sustain homogeneity across different regions. How frequent and continuous interactions could maintain the morphological variations revealed in this study should be explored in future work.

Second, our discussions have focused on PCs with a contribution rate of higher than 5 percentage points, which covers less than 40% of total cumulative proportions. Although the results of PCA on 3D data typically tend to be dispersed across more PCs [55,56,81], we should be careful not to overestimate the present results. The question of how frequently the Jōmon people interacted should be explored from additional aspects, including quantitative investigations of relevant archeological remains and mathematical simulations or modeling of the Jōmon population.

Third, some results from genetic or isotopic research on Jōmon skeletal remains, suggesting that Jōmon people interacted with populations from different regions, are consistent with our interpretation [82,83]. Cooke et al [84] argued that Jōmon people on the mainland of the Japanese archipelago were isolated from other islands, which might have contributed to the morphological homogeneity of Jōmon people. Nevertheless, to strengthen or verify this hypothesis, further genetic data and morphological data from other relevant regions overseas should be collected and combined [53,85,86]. For instance, our sample size is much larger than that in a study by Buck et al [53], which claimed that dietary factors influenced the shape of the Jōmon neurocranium, while our results should also be combined with isotopic dietary data in future work.

Fourth, while the dataset comprising larger sample sizes is acknowledged, potential sampling bias may exist, particularly due to restrictions on obtaining samples, notably in the Kyushu region. However, previous biometric studies have suggested that individuals from the northern Kyushu region do not exhibit significant differences compared to other regions [42]. This contributes to the overall reliability of the findings presented in this study. Our study exclusively focuses on cranial 3D data, although previous investigations have examined other parts of skeletal remains, such as limbs [36,39], warranting a comprehensive examination of their 3D data.

Finally, the selection of landmarks in this study is a critical aspect that warrants consideration. Although the chosen landmarks were selected for their clarity and suitability, the effectiveness of the landmarks in capturing morphologically significant changes may vary depending on various factors such as geographic regions, time periods, age groups, and other contextual factors.