Hominoid intraspecific cranial variation mirrors neutral genetic diversity.

Natural selection, developmental constraint, and plasticity have all been invoked as explanations for intraspecific cranial variation in humans and apes. However, global patterns of human cranial variation are congruent with patterns of genetic variation, demonstrating that population history has influenced cranial variation in humans. Here we show that this finding is not unique to Homo sapiens but is also broadly evident across extant ape species. Specifically, taxa that exhibit greater intraspecific cranial shape variation also exhibit greater genetic diversity at neutral autosomal loci. Thus, cranial shape variation within hominoid taxa reflects the population history of each species. Our results suggest that neutral evolutionary processes such as mutation, gene flow, and genetic drift have played an important role in generating cranial variation within species. These findings are consistent with previous work on human cranial morphology and improve our understanding of the evolutionary processes that generate intraspecific cranial shape diversity within hominoids. This work has implications for the analysis of selective and developmental pressures on the cranium and for interpreting shape variation in fossil hominin crania.

Synthesizing human cranial data and population genetic data has demonstrated that cranial variation tracks population history in a manner analogous to genetic data at a global scale (1–4). Most human cranial diversity is found within Africa, and diversity declines with increased distance from Africa. This follows patterns found in microsatellite and SNP data and is interpreted as evidence that nonselective evolutionary factors played a strong role in structuring the geographic pattern of human cranial diversity evident today (5–9). This finding is noteworthy because traditionally a greater emphasis was placed on natural selection and developmental pressures as factors shaping cranial variation in humans (10–12). Now, large population genetic datasets and quantitative genetics approaches have enabled more rigorous testing for both selective pressures and neutrality on different parts of the human cranium. The temporal bone and basicranium consistently reflect neutral patterns similar to microsatellite data (8, 13). However, the face shows evidence of climatic adaptation in certain human populations in which cranial distances exceed what would be expected under neutrality and instead correlate with temperature variables (2, 7).

Quantitative genetics provides a theoretical framework for testing for neutrality in phenotypic traits in a population (14, 15). This approach uses analyses derived from evolutionary genetics to detect departures from neutral patterns in phenotypic data. Globally, cranial diversity fits expectations of iterative founder effects, with similarity between populations decreasing exponentially as the geographic distance between them increases (16). Due to patterns of dispersal and expansion, human genetic diversity decreases with increasing distance from Africa. Distance from Africa also accounts for a portion of heritable variation in craniometric measurements (17). These patterns have been explored using both linear cranial measurements and 3D landmark data and with genetic models such as isolation by distance (IBD) and measures of population differentiation such as the fixation index (FST). Human postcranial elements have also been analyzed together with population genetic data, which shows pelvic variation follows geographic patterns with neutral genetic data. In contrast to the pelvis, long bone variation shows signatures of climatic adaptation (18–21). A quantitative genetics approach has also been applied to clarify the evolutionary forces driving hominoid cranial differences between species. For example, divergence time estimates for humans and Neanderthals that are calculated from cranial measurements are congruent with genetic divergence times and therefore imply that cranial form is evolving neutrally in these two lineages (22). However, extant hominoid cranial divergence has been characterized by stabilizing selection, which shows that patterns of cranial covariance structure and developmental integration are conserved across hominoids, with a few exceptions. These include cranial divergence in the lineages leading to Pongo, Hylobates, and Gorilla beringei beringei, which do not depart from neutral expectations derived from genetic divergence patterns (23).

Here, we look at cranial and genetic data together to determine the strength of the relationship between intraspecific cranial shape variation and neutral genetic diversity across 12 living ape species. Living apes are characterized by different amounts of both mitochondrial and nuclear genetic diversity because of their independent and varied population histories (24–26). These taxa also exhibit different magnitudes of neurocranial and facial shape diversity (27–29). If selective or developmental pressures are the primary factors influencing cranial variation within each taxon, we would not expect the magnitude of cranial shape variation to scale with neutral genetic diversity across all taxa. Stabilizing selection or developmental integration could constrain shape variation within a species instead of accumulating random differences through time together with neutral genetic loci. Conversely, environmental variables such as climatic differences and mechanical strain from mastication could inflate variation in certain taxa relative to others, thereby disrupting a relationship between neutral genetic variation and morphological variation across all taxa. Here, population genetic data from neutral autosomal loci reflect processes such as mutation, genetic drift, and fluctuations in population size through time and therefore provide a null hypothesis for testing for departures from these patterns in cranial data.

Results and Discussion

Using 396 adult crania from 12 hominoid taxa, we measured cranial variation with two different metrics, the average pairwise Procrustes distance (PPD) and the sum of eigenvalues (SEV) (). The first measure, PPD, reflects cranial shape differences among members of the same taxon, while SEV measures total shape variation within a taxon. Morphological data consisted of 34 homologous cranial landmarks divided into three units: (i) whole cranium, (ii) cranial vault (neurocranium), and (iii) face. Genetic data included 11 homologous noncoding autosomal loci across the same 12 taxa. Genetic variation was measured with pairwise nucleotide diversity (π), the number of segregating sites (θ), and the effective population size (N) ().

We performed 72 ordinary least squares regressions to determine the strength of the relationship between genetic and morphological diversity (). Twenty-four regressions were performed for each of the three landmark sets. Within each landmark set, the data were divided into a mixed-sex sample, females only, and males only. For each landmark set, within each sex category, eight regressions were performed to test all combinations of the different genetic and morphological metrics of variation: two with π (SEV, PPD), two with θ (SEV, PPD), two with N from π (SEV, PPD), and two with N from θ (SEV, PPD).

With a significance threshold of P value < 0.05, 57 of 72 regressions are statistically significant. If we apply a more conservative significance threshold of P value < 0.01, 42 of 72 regressions are significant. If we apply the Bonferroni correction to account for the multiple regressions that were performed (0.05/72) P value < 0.00069, 12 of 72 regressions are significant. Results reported within the text below use the P value < 0.01 significance threshold; results from all 72 regressions can be found in .

In the mixed-sex sample for the whole cranium, π accounted for 80% of the variance with SEV (P = 0.0001) (Fig. 1). In the female whole-cranium dataset, no results were significant at the 0.01 P value threshold. In the male whole-cranium dataset, π accounted for 77% of cranial variance with SEV (P = 0.0002).

Fig. 1.

Relationships between nucleotide diversity (π) and cranial shape variation (SEV) in the mixed-sex sample. (A) Whole cranium. (B) Cranial vault. (C) Face. Gg, G. gorilla gorilla; Hk, H. klossii; Hm, H. moloch; Hp, H. pileatus; Hsp, H. sapiens; Pa, P. abelii; Ppn, P. paniscus; Ppy, P. pygmaeus; Pts, P. troglodytes schweinfurthii; Ptt, P. troglodytes troglodytes; Ptv, P. troglodytes verus; Ss, S. syndactylus.

For the cranial vault, all regressions were statistically significant (P value < 0.01) for the mixed-sex sample and for males. For the mixed-sex cranial vault results, π accounted for 61% of the variance with SEV (P = 0.0028) (Fig. 1). For females, π accounted for only 47% of the variance with SEV (P = 0.0143). In the male cranial vault dataset π accounted for 61% of variance with SEV (P = 0.0026).

For the face, π accounted for 80% of the variance with SEV (P = 0.0001) in the mixed-sex sample (Fig. 1). None of the eight regressions of the female-only facial dataset were significant. For males, only two of eight regressions with the facial dataset were significant, π and θ with SEV.

The cranial vault shows a statistically significant relationship with genetic data in 20 of 24 regression analyses (P value < 0.01), while the face does so only in 6 of 24 regression analyses ().

The male cranial datasets were statistically significant in 18 of 24 regressions, and the female cranial datasets were statistically significant in only in 4 of 24 regressions.

Between the two metrics of morphological variation, SEV and PPD, we see mostly agreement in which regressions are statistically significant, with two exceptions. In the cranial vault, females across all genetic metrics are significant with PPD but are not significant with SEV. In the face mixed-sex sample, all results are significant with SEV but are not significant with PPD ().

We find no differences in which regressions are statistically significant between the two different measures of nucleotide variation (π, θ). The same combinations of variables (sex, landmark set, and morphological metric) are statistically significant for both N and the raw measures of genetic variation (π, θ), with only one exception: SEV in the male faces dataset which is significant with π and θ but is not significant with N.

Natural Selection or Neutrality in Hominoid Cranial Evolution?

A number of important implications arise from these results. First, our results suggest that the population processes that generate genetic variation at neutral loci explain a portion of the magnitude of cranial variation within each taxon, especially in the cranial vault. This result is consistent with patterns found in humans and indicates that neutral cranial evolution may not be unique to humans but rather may be part of a broader pattern in extant hominoids. In modern humans, the face tracks population history less closely than the rest of the cranium. The face may be more influenced by environmental pressures such as climate and bone remodeling due to masticatory strains. In humans, populations that live in extremely cold environments show departures from neutrality in aspects of nasal morphology, cranial breadth, and vault size and shape (2, 7, 30). In the primate cranium, mechanical strain has been shown to inflate variation in phenotypically plastic regions, especially the mandible and the face (12). Thus, the less consistent relationship between facial shape and genetic variation here may be partially driven by mechanical stress and phenotypic plasticity in hominoids. The cranial vault, however, is not subject to the same functional strains from mastication, is less variable than the face and mandible, and shows a statistically significant relationship with genetic data in the majority of the analyses here.

Despite the congruence between genetic and cranial data here, it is unlikely that only neutral evolutionary forces are acting to impact intraspecific cranial variation within hominoid taxa. Neutral processes such as mutation and genetic drift act in concert with developmental and selective pressures, and disentangling the differential effects of these processes on morphological and genetic evolution remains a long-standing challenge in evolutionary biology (31). It is important to note that the intraspecific focus of our study differs from recent work looking at how drift and selection impact variation and diversification between species (22, 23). It is possible that neutral population processes and levels of genetic variation explain a portion of cranial variation observed within the taxa included here but that directional or stabilizing selection were the dominant forces in driving or constraining diversification between these taxa. For example, Weaver and Stringer (22) show that between subspecies of Pan cranial differentiation is constrained relative to their genetic divergence. This suggests that cranial divergence between subspecies of Pan may be under stabilizing selection or that there is less variation available for genetic drift to act on because of developmental or genetic constraints. Here, subspecies of Pan fit the overall pattern across hominoids, that cranial shape variation scales with genetic diversity, but they consistently fall below the regression line, indicating that their cranial diversity is lower than expected given their genetic diversity (Fig. 1). Further analysis, which incorporates morphological variation within and between hominoids in a quantitative framework, could clarify rates of cranial versus genetic change at different taxonomic scales.

Sexual Dimorphism.

All the male-only regressions show a statistically significant relationship with genetic data for the whole cranium and cranial vault, but not for the face. In females, however, none of the regressions were significant for the whole cranium or the face, and only half were significant for the cranial vault. In total, 18 of 24 regressions were significant in males, and only 4 of 24 were significant in females (). This result is interesting in the context of other work that demonstrates that morphology reflects phylogeny more in males than females, especially in highly dimorphic taxa (32). Here, males show a more consistent relationship with genetic data than do females due to the differences in male and female cranial variation within the same taxon. For example, in highly dimorphic taxa (Pongo and Gorilla), females show less cranial shape variation than males, but this is not the case for less dimorphic taxa such as Homo sapiens and Pan paniscus (Fig. 2). The stronger relationship between genetic and cranial data in males, coupled with females showing less cranial shape variation than males, suggests that selection may be playing a role in constraining variation in female cranial shape relative to males in certain taxa.

Fig. 2.

Cranial shape variation within extant apes (SEV) in the mixed-sex sample and in females and males separately. Taxa are listed in descending order by whole-cranium mixed-sex sample. (A) Whole cranium. (B) Cranial vault. (C) Face. Abbreviations as in Fig. 1.

Additionally, Gorilla and Pongo are not only the most sexually dimorphic apes; they also show some of the highest levels of genetic diversity. Therefore, it is key to note the single-sex results in these taxa particularly. If we look only at the single-sex whole-cranium analyses (SEV), Pongo and Gorilla are highly variable compared with other taxa, especially in the male-only datasets. These species also show the highest π values of all of the apes. If cranial dimorphism were inflating cranial variation in the mixed-sex sample and driving the relationship between genetic and cranial variation in these taxa, we would not expect single-sex data to show the statistically significant relationship with genetic data that they do. Additionally, the mixed-sex gorilla data (SEV) show more cranial shape variation than would be expected, given their genetic diversity; this is due to marked sexual dimorphism in this species, especially in the cranial vault (Fig. 1).

Hominoid Population History.

If population history has impacted both genetic and morphological intraspecific variation in apes, the central question then becomes: What demographic and ecological factors drove this parallel change in molecular and skeletal diversity? Fluctuations in population size and structure through time, which are gleaned from genetic data, together with biogeographic information, provide baseline explanations for our findings in cranial data. All living apes have smaller population sizes than humans, but genetic diversity has been maintained in many species as a vestige of large ancestral population sizes, population substructuring, and older lineage ages (25). The reverse is true in modern humans, in which a recent origin and a population size reduction followed by rapid expansion has resulted in lower genetic diversity (Table 1) (33, 34).

Table 1.

Cranial shape variation: SEV and PPD in mixed-sex, whole-cranium datasets, pairwise nucleotide diversity (π), geographic range, and population history for extant hominoids

Hominoid taxon	N	n	SEV	PPD	π (%)	Census size	Geographic distribution	Population history inference
P. abelii	8 f, 10 m	6	0.0200	0.0166	0.42	∼7,300	Sumatra	Long-term fragmented range, recent range reduction (25, 40)
P. pygmaeus	20 f, 15 m	10	0.0142	0.0091	0.35	∼50,000	Borneo	Long-term fragmented range, recent range reduction (25, 40, 41)
G. gorilla gorilla	29 f, 41 m	14	0.0141	0.0102	0.19	∼95,000	Central Africa	Constant population size, almost continuous distribution, recent range reduction (25, 38, 39)
S. syndactylus	17 f, 22 m	6	0.0090	0.0103	0.21	∼190,000	Sumatra, Malay Peninsula	Long-term widespread, shrinking but continuous populations (43)
H. sapiens	18 f, 20 m	100	0.0065	0.0098	0.10	∼7 billion	Cosmopolitan	Bottleneck followed by recent massive range expansion (33, 34)
P. troglodytes schweinfurthii	3 f, 8 m	10	0.0057	0.0062	0.15	∼89,000	Congo River to W. Uganda, Rwanda, W. Tanzania	Bottleneck, expansion, and recent range reduction (25, 35, 36)
P. troglodytes troglodytes	50 f, 26 m	10	0.0056	0.0064	0.18	∼90,000	Central Africa Sanaga River to Congo River	Constant population size, recent range reduction (25, 35, 36)
H. moloch	6 f, 8 m	4	0.0050	0.0080	0.21	∼2,500	Java	Constant island population, recent range reduction (43)
P. troglodytes verus	12 f, 10 m	10	0.0044	0.0064	0.12	∼55,000	West Africa, Senegal to Nigeria	Bottleneck, expansion and recent range reduction (25, 35, 36)
P. paniscus	21 f, 17 m	9	0.0037	0.0057	0.09	∼50,000	Central Africa, South of Congo River	Bottleneck and continuous restricted range (25, 36)
H. klossii	10 f, 11 m	2	0.0035	0.0062	0.08	∼25,000	Mentawai Islands	Long-isolated island populations (44)
H. pileatus	3 f, 3 m	4	0.0031	0.0052	0.05	∼40,000	SE Thailand, Cambodia, SW Laos	Restricted range (44)

Taxa are ordered by SEV. N = number of cranial samples in females (f) and males (m); n = number of genetic samples. All population census size estimates are from the International Union for the Conservation of Nature Red List.

Within Pan, π and N are highest in Central chimpanzees (Pan troglodytes troglodytes), intermediate in Eastern chimpanzees (Pan troglodytes schweinfurthii), and lowest in Western chimpanzees (Pan troglodytes verus) and bonobos (P. paniscus) (35). Data from Bayesian population modeling in chimpanzees suggests that Eastern and Western populations experienced a bottleneck after their divergence before expanding to their current range. In contrast, Central chimpanzees show a recent range expansion without evidence of a bottleneck (36). Bonobos and common chimpanzees were separated by the formation of the Congo River ∼1.5–2 Ma. This barrier inhibited gene flow and restricted bonobos to a small area south of the river. Periodic contractions of forest cover in this region may have forced bonobos into a bottleneck, which is consistent with their low genetic diversity and low cranial shape diversity (37). Nucleotide diversity within Western lowland gorillas (G. gorilla gorilla) is close to estimates within Central chimpanzees and is higher than in other members of the genus Pan (Fig. 3). During the Last Glacial Maximum, African rainforests became fragmented, a process that was reversed postglacially when forest patches expanded and rejoined (38, 39). Gorilla populations may have become discontinuous during this time, creating separate reservoirs of diversity. This population structure would have resulted in the maintenance of genetic (and phenotypic) diversity by providing novel mutations when groups resumed gene flow.

Fig. 3.

Nucleotide variation (π) in hominoids. Average π values are shown for each locus (). Abbreviations as in Fig. 1.

Orangutans show the highest levels of nucleotide diversity among the great apes, with the Sumatran species (Pongo abelii) being more variable than the Bornean species (Pongo pygmaeus) (40, 41). Sumatran orangutans have three deeply structured genetic clusters, indicating long-term separation of these groups (42). Our results support a complex population history for orangutans that is marked by high intraspecific morphological and genetic diversity among hominoids despite small population sizes. In comparison with orangutans, hylobatids are more species rich and geographically continuous. There is evidence of a recent radiation of hylobatid species less than 2 Ma, followed by continued gene flow between certain species (43, 44). These processes have reduced variation between species but may have acted to maintain variation within certain species. Here, nucleotide diversity is higher within Symphalangus syndactylus and Hylobates moloch than within members of the genus Pan. S. syndactylus also shows greater cranial shape diversity than Pan. S. syndactylus still maintains a geographic distribution on Sumatra and a central section of the Malay Peninsula and has a large census size relative to other hylobatids. Genetic variation in H. moloch is higher than might be expected given its critically low census size of 2,500 individuals, although its current distribution in forest fragments in Western and Central Java is not representative of the historical range of the species. Hylobates klossii lives exclusively on the Mentawai islands and has the smallest geographic range of the hylobatids in this study; its genetic and cranial and variation are among the lowest presented here. Hylobates pileatus, in southeastern Thailand, Cambodia, and southwestern Laos, shows the lowest genetic diversity estimates and cranial variation of all hominoid species here.

Implications for Fossil Hominins.

The finding that cranial morphology preserves signals of past population history in hominoids can also guide our understanding of variation in fossil hominin crania. For example, extant hominoids often serve as modern analogs of variation which inform inferences of intraspecific variation in extinct groups. Our results suggest that the population history of extant apes should be considered when choosing analogs of intraspecific variation for fossil hominins. Accordingly, modern H. sapiens may be a suboptimal analog for variation in fossil hominins, despite their close phylogenetic relatedness. Modern humans have a unique population history featuring at least one severe bottleneck followed by rapid expansion and repeated founder events (33, 34). If these features of human population history have a major impact on cranial variation, then modern human variation provides a limited model of cranial variation for extinct hominins. For example, for fossil species such as Homo erectus—a cosmopolitan hominin species with a temporally longer lineage than modern humans—we might expect more variation than we see in humans today if this species did not experience an equivalent bottleneck and rapid expansion. Additionally, for fossil hominin populations from the same species and time horizon (e.g., Homo naledi from the Dinaledi chamber, South Africa, H. erectus from Dmanisi, Georgia, and the fossil hominins at Sima de los Huesos, Spain), the amount of cranial variation in adults may broadly reflect population genetic structure (45–47). Overall, results here suggest that intraspecific morphological variation in hominin crania can be viewed through a population genetics framework—with consideration of how multiple different population models could explain the observed levels of variation.

Future Research.

The results we report here open possibilities for future analyses within the primates, including those at different geographic sampling scales and with the inclusion of additional skeletal elements such as postcrania and dentition. This work also may have implications for developmental analyses of the cranium. For example, if neutral genetic diversity explains a portion of intraspecific cranial shape variation, then studies assessing variation in different developmental modules of the cranium could account for different population histories when comparing taxa.

More broadly, this work provides preliminary empirical support that neutral population processes have impacted extant hominoid cranial morphology and evolution, and that this pattern may be relevant to other taxonomic groups.

Materials and Methods

Morphological Data.

A total of 396 adult crania from 12 taxa (species and subspecies) were included in this analysis: H. sapiens, P. paniscus, P. troglodytes troglodytes, P. troglodytes verus, P. troglodytes schweinfurthii, G. gorilla gorilla, P. pygmaeus, P. abelii, S. syndactylus, H. moloch, H. klossii, and H. pileatus (Dataset S1).

Morphological data consisted of 34 homologous cranial landmarks (each consisting of a set of x, y, z coordinates) capturing cranial shape differences () (48, 49). All landmark data were subjected to a generalized Procrustes Analysis (GPA) to project them into a common shape space. The GPA superimposes the centroids of each individual’s landmark configuration, then scales all landmark configurations to unit size and rotates all specimens around that centroid. Differences in translation, size, and orientation are eliminated during this step, so that only differences in shape remain among specimens. Landmark data were divided into three analytic units: (i) whole cranium, (ii) cranial vault, and (iii) face, consisting of 34, 12, and 22 landmarks, respectively. Three separate GPAs were performed for each of the three landmark units.

We applied two types of morphological analyses: average PPD and SEV. The average PPD for each taxon was calculated as the mean squared distance between all pairs of individuals belonging to the same taxon optimally aligned in Kendall’s shape space () (50, 51). The SEV for each taxon was calculated from the variance–covariance matrix of the superimposed coordinates calculated separately for each taxon. This is a symmetric square matrix in which the diagonal elements are the variances of the individual shape coordinates and the off-diagonal elements are the covariances among coordinates subsequent to the GPA superimposition described above. This value equals the cumulative variance in that group across all landmarks (Fig. 2 and ) and is also equivalent to the group’s Procrustes variance, which measures the mean squared Procrustes distance of each specimen to the average shape (52). Quantitative genetic theory predicts a linear relationship between SEV and neutrally evolving population genetic data; however, this is not the case for PPD (53). We chose to include both measures of morphological variation here to provide evidence that the patterns of cranial shape variation within taxa are similarly robust when different methods are used.

Tests were performed to determine how sample size influenced the mean PPD within each group. The largest sample size was 76 individuals for P. troglodytes troglodytes. Two individuals from the same species were randomly sampled, and the PPD between them was calculated. The resampling procedure was repeated 10,000 times with replacement for each taxon, and the average pairwise distance was recorded for random subsets of 75, 50, 20, 10, and 5 P. troglodytes troglodytes individuals. For each different sample size, the resampling procedure was repeated 10 separate times and then was averaged. Across the different sample sizes, the average PPD values were 0.00641 for all three of the largest sample sizes (n = 75, n = 50, and n = 20), 0.00645 for n = 10, and 0.00642 for n = 5. Varying the sample sizes did not yield appreciably different average PPD values, thus confirming that sample sizes used here were adequate for capturing intraspecific cranial shape diversity that reflects a larger taxon-wide trend. This is especially relevant for estimation of PPD for single-sex, within-taxon samples, which were represented by the smallest number of individuals.

Sample size tests for SEV were also performed. Overall, the smallest mixed-sex cranial dataset was six individuals (H. pileatus). All taxa were sampled down to six individuals (three females, three males) for the whole-cranium dataset, and SEV was calculated (). Even with sample sizes of six individuals, we see similar trends in the magnitude of variation within taxa. Pongo, Gorilla, and Symphalangus remain the most variable taxa, and P. paniscus, H. klossii, and H. pileatus are the least variable. We see the largest difference in SEV values between the sample size of six and the larger sample size of 76 in P. troglodytes troglodytes. This taxon represents the largest sample size in our analysis.

Genetic Data.

Genetic diversity was summarized by pairwise nucleotide diversity (π) and number of segregating sites (θ) using the software SITES (54). Sequence data from 11 homologous noncoding autosomal loci across all 12 taxa were downloaded from GenBank (). The genetic estimator of π is calculated by randomly sampling two individuals within each taxon, with replacement, and then taking the average nucleotide differences between pairs. This was performed for each locus separately and then averaged across all loci to arrive at a single π value for each taxon (Fig. 3 and ). Nucleotide diversity estimates (π, θ) from the autosomal loci chosen here reflect neutral evolutionary processes and serve as a proxy for genome-wide impacts of population history.

Effective population size was calculated using the standard equations: N = π/(4μ) and N = θw/(4μ). The mutation rate (μ) was derived from the average pairwise differences between two species (for all loci) divided by the number of generations between the estimated divergence of the two species () (55). The following generation times were used in the N calculations: H. sapiens, 29 y; Pan, 25 y; Gorilla, 19.3 y; Pongo, 26.7 y; and Hylobates and Symphalangus, 20 y. The following divergence times were used: Homo–Pan, 7 Ma; Pan–Gorilla, 13.5 Ma; Gorilla–Pongo, 16 Ma; Pongo–Symphalangus, 20 Ma; and Pongo–Hylobates, 20 Ma (55).

Regression Analyses.

Finally, 72 ordinary least squares regressions were performed to determine the strength of the relationship between genetic and morphological diversity (). Ordinary least squares (OLS) regressions were performed with raw π and θ values as well as N and N. Each of these four genetic values was separately regressed against the corresponding SEV and PPD values. This was done for each landmark set and for both sexes together and separately. Additionally, phylogenetic generalized least squares (PGLS) regressions were performed to determine if the results from the OLS regressions were the result of close evolutionary relationships between taxa. A tree file of all species used in this analysis was generated from the website https://10ktrees.nunn-lab.org/index.html. This file was loaded into R, and PGLS regressions were performed using the packages (caper) and (ape). All results showed that the relationship between genetic and cranial variation was not a result of close evolutionary relationships between the taxa sampled here, with lambda values at or close to 0 (indicating that the data are not compatible with a Brownian motion model of evolution).