Geometric Morphometrics Provides an Alternative Approach for Interpreting the Affinity of Fossil Lizard Jaws

The study of fossils provides a unique window for research into the evolutionary history of taxa. The age, geographic origin, and palaeoenvironment of fossils provide important sources of evidence for the evolution of morphological characters, past distribution of taxa, and wider environmental changes. This is particularly true for recent fossils (<500 thousand years [kyr]) that may have close living relatives, or even represent living species, with ecological tolerances that are well understood. For such fossils to be of any use, however, they first must be identified to some taxonomic level with confidence. The reliability of the alpha taxonomy is of high importance because it often is used in broader studies to quantify past changes in diversity or to constrain molecular divergence analyses (e.g., Parham et al., 2012; Bell et al., 2010; Mannion et al., 2015; Slavenko et al., 2016).

The Holocene-Pleistocene fossil reptile assemblages of Australia have not been well studied in comparison with their mammal counterparts (e.g., Travouillon et al., 2006). Members of reptile clades can be found in samples from many fossil localities, potentially representing most of the major components of an exceptionally diverse living squamate fauna; however, interpretation is inhibited by a poor understanding of reptile osteology at low taxonomic levels and an inability to make objective comparisons. Variation within and between species tends to involve subtle differences in the shapes of processes and in relative proportions that can be a challenge to compare holistically (Evans, 2008; Hollenshead et al., 2010; Sherratt et al., 2015). Specific characters can be defined and used as apomorphies (e.g., Hutchinson, 1997), but this requires some baseline knowledge of variation within the taxon being examined, and such data may not be available (Bell and Mead, 2014). Many skeletal collections lack adequate samples of lizards, and the published descriptions and images can be of limited use because they are most-often focused on the articulated cranium as a whole. The taxonomy of modern species tends to be derived from analysis of genetic data and external characters (e.g., scale number, proportions) and does not tend to provide any information on osteological characters (e.g., McLean et al., 2013). Authors documenting fossils tend to provide outline drawings or photographs of specimens (e.g., Covacevich et al., 1990; Hocknull et al., 2007; Prasad and Bajpai, 2008) which assist in broad comparisons but cannot convey the full three-dimensional aspect of the bones. Quantitative measurements representing tooth shape can provide more-objective comparisons (e.g., Jones et al., 2009; Hollenshead et al., 2010), but our understanding of the extent of associated variation in all taxa is poorly developed. Moreover, it may not be possible to measure the teeth accurately in a fossil because of tooth wear or breakage.

Geometric morphometrics provides an alternative and potentially more-objective approach for characterizing and comparing bones (Adams et al., 2004), although it has yet to be widely applied to the isolated bones of Quaternary reptile fossils. Comparing anatomical structures generally involves the use of landmarks and outlines to quantify variation amongst specimens, and analyses require multivariate statistics. Geometric morphometrics is considered superior to previous forms of biometrics (e.g., isolated linear measurements) because it records the geometric relationship between a cloud of landmarks and/or set of curves, providing a holistic measurement of overall shape in contrast to isolated linear measurements (Adams and Collyer, 2009). This facilitates standardization between specimens of different size via Procrustes superimposition and permits visualizations that aid interpretation. Over the past decade, geometric morphometrics has become an increasingly accessible approach for morphological analyses and has been used to characterize and compare shape variation amongst two-dimensional images of reptile skulls (e.g., Stayton, 2005; Jones, 2008; Meloro and Jones, 2012; Sanger et al., 2013; Fabre et al., 2014; Openshaw et al., 2016); however, it also has been applied to three-dimensional reptile anatomy using X-ray computed tomography (CT) (Parr et al., 2012; McCurry et al., 2015). In at least one case, it has been used as an approach to associate fossils with modern taxa (Dollion et al., 2015).

Here we use three-dimensional geometric morphometrics (see Zelditch et al., 2012) to more-objectively characterize the maxillae of modern agamid lizards from southern Australia to facilitate comparisons with a recent (<20 kyr) fossil specimen from Kangaroo Island. Ctenophorus decresii is the only agamid species present on Kangaroo Island today, but fossil remains provide the opportunity to examine whether different taxa were present there in the past.

We analyzed an almost-complete agamid fossil maxillae (P53917) from Kelly Hill Caves on Kangaroo Island (McDowell et al., 2013). The specimen was recovered from sediments with an age range of 11,645–10,360 yr before present (95% confidence interval) according to a chronological model developed by Matthew C. McDowell using Bayesian analysis. The fossil assemblage was carefully excavated layer by layer and dated using U-Th dating of speleothems, accelerator mass spectrometry (AMS) radiocarbon-dating of bone, and optically stimulated luminescence of quartz grains. Excavated sediment was wet-sieved using 1.5-mm mesh. The residues of small vertebrate remains were dried and then sorted (picked) for taxonomically identifiable specimens. The specimen examined here is entirely removed from the matrix and has excellent surface preservation, clearly showing the location of foramina and sutural facets (Fig. 1). The specimen is essentially complete except that the distal edges of the posterodorsal process may be rounded and a portion of the dorsal edge of the facial process is broken and missing.

We characterized 10 modern species, each represented by a single specimen: Amphibolurus muricatus (Museum of South Australia [SAMA], R21375), Amphibolurus norrisi (SAMA R60767), Ctenophorus decresii (SAMA R28618), Ctenophorus fordi (SAMA R34489), Ctenophorus pictus (SAMA R28608), Pogona barbata (SAMA R32503), Pogona vitticeps (SAMA R18545), Tymanocryptis lineata (SAMA R59721), and Rankinia diemensis (SAMA R269B). These taxa represent agamids living in South Australia today or that have been reported from South Australian fossil deposits (Owens and Graham, 2009; Reed and Bourne, 2009). Our set of specimens was intended to be a minimum sampling to assess the ease with which different species could be characterized and was used as a pilot data set for comparison of the modern fauna with recent fossils.

The fossil specimen was subjected to micro X-ray CT at Adelaide Microscopy using a Bruker Skyscan 1076 at a resolution of 9 μm. The CT scan reconstructions were obtained for each comparison species from the scanned specimen database at South Australian Museum. All CT scans used in this study were reconstructed using NRecon software provided by the scanner manufacturer (NRecon, version 1.6.9.4, Skyscan, Kontich, Belgium). The reconstructed scans were digitally segmented to extract the right maxilla from the cranium, and surface files were created using the “segmentation editor” in Avizo 8.1 software (VSG Inc., Burlington, Massachusetts, U.S.A.) with minimal smoothing (Fig. 1).

Landmarks were placed on the surface files of the fossil and comparative maxillae using the software package “IDAV Landmark” (Wiley et al., 2005). Twenty-two fixed landmarks and 30 semi-landmarks were placed at equivalent points on the three-dimensional surface files: i.e., at the most-extreme points of particular features and along the major curves (shown in Fig. 2). Each fixed landmark consisted of a single point on the surface of the bone with x, y, z coordinates, and each semi-landmark was part of a collection of 10 evenly spaced points along a curve. For the fossil specimen, missing landmarks (2 and 6) were estimated using a multivariate regression method (where each landmark with missing values is regressed on all other landmarks for the set of complete specimens, and the missing landmark values are predicted by this linear regression model) using the R package ‘geomorph' (Adams and Otarola-Castillo, 2013) in Rstudio (www.rstudio.com/products/rstudio/download2/).

We repeated the following analysis on two data sets: one data set contained all the extant comparison specimens and the other contained all the extant comparison specimens and fossil specimen P53917. Both the Procrustes superimposition and principal component analysis were carried out using the R package ‘geomorph'.

We performed a Procrustes superimposition to scale the sets of landmarks and fit the shapes to each other to remove the shape differences related to absolute size (but not shape allometry) (Rohlf and Slice, 1990). This effectively fitted the landmark constellations around a mean constellation. Two distance matrices containing the Procrustes distances between each specimen for both data sets were calculated using ‘geomorph' in R, which provided an explicit and quantitative measure of overall shape similarity. Semi-landmark tangent sliding directions were specified using the “Procrustes distance criterion” (Bookstein, 1997). Next, we performed a principal components analysis (PCA) to find the linear combination of variables that represent maximum variance within the current specimen sample. A Mantel test (Mantel, 1967) was performed in R to evaluate the similarity of the Procrustes distances of the extant specimens resulting from analysis of both data sets.

We used ‘geomorph' to first produce a mesh whose shape was defined by the mean shape of its landmark coordinates. We then warped this mesh into the shapes represented by the maximum and minimum values of PC1 and PC2. This provided a visualization of the PCA results from which we could determine major shape differences for PC1 and PC2. All mesh warping was performed using the thin-plate spline method (Bookstein, 1989).

The Procrustes distances, shown in Table 1, indicate that the species pairs closest to one another in overall shape are A. muricatus–A. norrisi, A. muricatus–C. pictus, A. norrisi–P. barbata, and C. pictus–P. vitticeps. The two taxa most dissimilar to others are C. fordi and Tympanocryptis lineata, with relatively large Procrustes distances between them and other specimens included in the analysis. The species pairs most dissimilar to one another were T. lineata–A. norrisi and T. lineata–C. fordi. Ctenophorus decresii is most similar to C. pictus.

The PCA (Fig. 3) revealed that PC1 (34.9% of total shape variation) describes shape differences involving the anterior process, where it ranges from being relatively deep at negative values (e.g., A. norrisi) to relatively shallow at positive values (e.g., T. lineata). It also describes the shape of the facial process, from a swept-back process that is relatively broad at its base at more negative values (e.g., P. vitticeps) to an upright process that is narrow at its base at more positive values (e.g., T. lineata). The overall shape also contrasts between a dorsoventrally deep maxilla at negative values to a dorsoventrally shallow maxilla at more positive values. Maxillary shape associated with minimum and maximum values of PC1 are shown in Figure 4.

The PC2 (21.0% of variance) describes differences in shape of the suborbital margin, which is quite smooth at positive values (see A. norrisi), and more irregular at negative values (see C. fordi). The subnarial margin also contrasts considerably in shape, from being small and having an upright anterior margin at negative values to being large with a sloped anterior margin at positive values. The anterior process graduates from having a subtle bifurcation at negative values to having a pronounced bifurcation at positive values. The shape of the maxillae for this PC differs from being curved along the long axis at negative values to being relatively straight at positive values. Maxillary shape associated with minimum and maximum values of PC2 are shown in Figure 4.

The PC1 recognizes the difference mainly between T. lineata and other agamids whereas PC2 recognizes the difference mainly between species of Ctenophorus. Amphibolurus muricatus and A. norrisi have similar scores for PC1 and PC2. This also is true for P. barbata and P. vitticeps. The three species of Ctenophorus have similar scores for PC1. As indicated by the Procrustes distances, T. lineata and C. fordi are particularly different to everything else in the sample, plotting at extremes of PC1 and PC2, respectively. The PC3 (17.0%) recognizes the difference mainly between P. vitticeps and P. barbata whereas PC4 (9.7%) recognizes the differences between R. diemensis and all other agamids.

The Procrustes distances of the fossil specimen (Table 2) show that it is most similar in shape to A. muricatus and least similar to T. lineata. The fossil also is quite dissimilar to C. decresii (Table 2). A Mantel test comparing the Procrustes distances between extant specimens for both data sets produced an observed correlation coefficient of 0.999 and a P-value of 0.001 (after 999 permutations). These results indicate that inclusion of the fossil specimen does not greatly alter the distances between the modern specimens. The pairs most similar (A. muricatus–A. norrisi) or dissimilar (T. lineata–A. norrisi) amongst the extant samples remained the same.

The fossil had scores for both PC1 and PC2 (together 51.6 % of total shape variation, Fig. 5) that were similar to those of the two species of Amphibolurus having a relatively deep anterior process, swept-back facial process that was wider ventrally, and deeper in shape dorso-ventrally. With respect to PC3 (15.6%), the fossil had similar scores to C. pictus and C. fordi and, for PC4 (9.6%), it had similar scores to R. diemensis. Inclusion of the fossil in the PCA had little effect on the shape differences described by the main axes of shape variation (PCs 1–4) or the distribution of the extant species on those PCs. The fossil did not plot near C. decresii on any of the four main axes of variation (76.8%).

The fossil maxilla does not have a shape consistent with that of C. decresii, the only agamid currently present on Kangaroo Island. Instead, it most closely resembles the maxilla of species of Amphibolurus. One species of this genus, A. norrisi, currently occurs on two mainland areas immediately adjacent to Kangaroo Island (southern Yorke and Eyre Peninsulas), suggesting a good habitat match for this species with Kangaroo Island. If this was the species present, its loss from the island could be related to isolation of the island at 9,000 kya (Adams et al., 2016), possibly because of the cooling effect of the surrounding ocean or as a consequence of the island area effect (MacArthur and Wilson, 1967). Confident allocation of the fossil specimen will require greater sampling of the living candidate species (A. norrisi and A. muricatus) both to ensure that these two can be distinguished from each other and to demonstrate to a higher degree of probability that the fossil falls within their range of variation. Similarly, a larger sample of C. decresii would allow us to better test that its range of variation would not accommodate the fossil. We note, however, that even single individuals of the candidate species appear to provide a useful framework for comparing the shape of the fossil and evaluating its likely identity. Our results provide the first evidence that in the past 20,000 yr, Kangaroo Island supported a different and possibly more-diverse reptile fauna.

Geometric morphometric analysis of three-dimensional scans provides a very promising approach for more-objective comparisons of microvertebrate fossils. Previous geometric morphometric analyses of reptiles were carried out on the complete cranial structure, but we show here that a similar approach may be used for the analyses of the isolated elements often recovered in fossil deposits. Although the sample size in our pilot study limits our interpretations, we still are able to provide a repeatable, explicit, and quantitative index of the overall shape differences between specimens and can examine the relative distribution of the specimens in an empirical shape space. A larger sample size likely could allow much greater discrimination power (Close and Rayfield, 2012; Wilson et al., 2013). This approach would be a potentially powerful tool for analysis of new fossils and evaluation of previous identifications (e.g., Holman and Case, 1988; Hsiou et al., 2016).

The method remains somewhat subjective with respect to landmark choice and sample used but, once established, it allows explicit comparisons in terms of similarities and differences in shape. We chose landmarks with the aim of providing the best and most balanced representation of the bone shape and also so that they may be placed unambiguously on every specimen under the same definition (see Appendix 1). Although we made our best effort to get the most-likely species for comparison, we still may have missed some taxa contained in the deposits. Nevertheless, the extreme unlikelihood of some matches may still be indicated by the analysis. Although we acknowledge that a simple equivalence between a fossil and a similarly shaped living species risks overinterpretation (Bell et al., 2010), nevertheless, the most-successful comparisons likely will be from younger sites where the number of extinct species probably is lower. Hence, geometric morphometrics adds another tool to reduce the subjectivity and “covert biases” noted by Bell et al. (2010).

The ability to estimate missing landmarks removes a long-standing restriction of geometric morphometrics (Gunz et al., 2009; Arbour and Brown, 2014). Before that development, missing landmarks needed to be removed from all specimens in the analysis, effectively excluding much of the shape information captured by the respective landmarks. This approach maximizes the morphometry of all specimens, increasing the power to discriminate between species despite within-species variance.

We thank R. Williams and A. Parker Watson (both University of Adelaide) for their expertise in computed tomography. MJ is currently supported by Australian Research Council Discovery Early Career Researcher Award DE130101567.

Images of a surface model of an example maxilla shown in labial view (A), lingual view (B), anterior view (C), to show components included in landmark definitions, and table containing the definitions of all landmarks used in the analyses.