Skip to main content

Phylogenetic conservatism in skulls and evolutionary lability in limbs – morphological evolution across an ancient frog radiation is shaped by diet, locomotion and burrowing



Quantifying morphological diversity across taxa can provide valuable insight into evolutionary processes, yet its complexities can make it difficult to identify appropriate units for evaluation. One of the challenges in this field is identifying the processes that drive morphological evolution, especially when accounting for shape diversification across multiple structures. Differential levels of co-varying phenotypic diversification can conceal selective pressures on traits due to morphological integration or modular shape evolution of different structures, where morphological evolution of different modules is explained either by co-variation between them or by independent evolution, respectively.


Here we used a 3D geometric morphometric approach with x-ray micro CT scan data of the skull and bones of forelimbs and hindlimbs of representative species from all 21 genera of the ancient Australo-Papuan myobatrachid frogs and analysed their shape both as a set of distinct modules and as a multi-modular integrative structure. We then tested three main questions: (i) are evolutionary patterns and the amount and direction of morphological changes similar in different structures and subfamilies?, (ii) do skulls and limbs show different levels of integration?, and (iii) is morphological diversity of skulls and limbs shaped by diet, locomotion, burrowing behavior, and ecology?.


Our results in both skulls and limbs support a complex evolutionary pattern typical of an adaptive radiation with an early burst of phenotypic variation followed by slower rates of morphological change. Skull shape diversity was phylogenetically conserved and correlated with diet whereas limb shape was more labile and associated with diet, locomotion, and burrowing behaviour. Morphological changes between different limb bones were highly correlated, depicting high morphological integration. In contrast, overall limb and skull shape displayed semi-independence in morphological evolution, indicating modularity.


Our results illustrate how morphological diversification in animal clades can follow complex processes, entailing selective pressures from the environment as well as multiple trait covariance with varying degrees of independence across different structures. We suggest that accurately quantifying shape diversity across multiple structures is crucial in order to understand complex evolutionary processes.


Understanding morphological evolution, and the underlying mechanisms that generate the enormous phenotypic diversity we see, is a central aim in evolutionary biology [1,2,3,4]. Phenotypic diversity often is correlated with ecology and behaviour, especially in traits for which form and function are tightly associated due to evolutionary and ecological pressures [5,6,7,8]. However, while some clades display extensive ecological and morphological variation that is correlated with lifestyle, others retain ancestral environmental niches and conserved body shape patterns that are better explained by phylogenetic conservatism on a shared ancestral lifestyle [9, 10]. These differing patterns of diversification are best illustrated in related groups of species where one group might display more phenotypic diversification than another due to different selective pressures [4, 11]. There are many examples of this in the species-rich radiations of characiform fishes [11], gobies and cardinal fishes [12], passerine birds [13], archosaurs [14], and many others.

While diverse evolutionary processes can generate phenotypic change, morphological evolution is typically inferred from integration or co-variation among multiple traits [15]. Body shape patterns can usually be broken down into ‘modules’, which are characterized by more internal integration within them, than externally among them [15]. Therefore, each module displays a certain amount of independence from other modules and can differ developmentally, genetically, and in the way they respond to selection [15, 16]. While many phenotypic changes across a radiation are modular in this way [17], shape diversification can follow a more complex pattern of integrative co-variation between modules and show correlated morphological variation among them [15, 18]. The degree of shape-co-variation between modules is due to the interplay between morphological integration and modularity, where morphological modules evolve in concert with others and in which morphology evolves independently among different structures, respectively [15]. High morphological diversity could be correlated with modularity, as autonomy among different structural units might promote higher independent morphological changes due to the evolutionary lability necessary for adaptive shifts [19,20,21]. Conversely, morphological integration could be one of the causes leading to convergence among unrelated clades [22, 23]. Integration and co-variation among modules should also shape the morphological evolution of individual organisms, as some modules might be subject to strong selective pressures from the environment, whereas others might be phylogenetically constrained. Therefore, identifying the patterns of variation in each module, while accounting for integration among them, is crucial in order to study morphological evolution and the processes that might have driven it [23].

Due to the close relationship between form and function, some morphological traits are likely to be more closely linked to the ecology of an organism than others [24]. For example, Zaaf & Van Damme [25] proposed the idea of evaluating morphological differences between and within distinct modules in limbs, in relation to functional traits like locomotion, and tested it in climbing and ground-dwelling geckos. Limb shape might provide the most insight into the ecotype a species occupies, as it is closely correlated with its performance, and thus, locomotion through the environment [26, 27]. Similarly, Cornette et al. [28] looked at both the skull and mandible in shrews in order to disentangle the relationship between diet, ecological factors, and head shape evolution. On the other hand, some modules might be correlated with life history traits or not be under selection as functional traits [29]. Moreover, inferring adaptive processes by looking at the ‘wrong’ structure might be uninformative, and in some cases even misleading. Assessing morphological evolution in a group of organisms provides more valuable information when looking at a wide range of phenotypic traits, but may also increase the difficulty of data interpretation, due to complex co-variation processes between different structures.

Anuran amphibians are an ideal model group in which to investigate morphological evolutionary patterns: they display a highly derived morphology compared to other terrestrial vertebrates [30], yet their body plan has been relatively conserved since the early Jurassic [31, 32]. Despite phylogenetic constraints on their appendicular skeleton as an adaptation to saltatory locomotion [33], substantially different body shape patterns have evolved independently across several clades [34]. Frogs and toads have adapted to a wide array of extreme environments through a combination of behavioural, physiological, and morphological mechanisms. Extreme morphological shifts are usually associated with unique locomotor types, such as gliding in “flying” frogs [35], or with specialised locomotion, such as the improved swimming ability in frogs like pipids [36]. Similarly, strong shape changes are observed in burrowing frogs and toads that have adapted to desiccating conditions in arid and semi-arid environments [26, 37]. Morphological convergence in burrowing frogs has been documented across numerous clades, in both forward (head and forelimbs first), and backward (hindlimbs first) burrowing species, with backward burrowing being the most common digging type in frogs and toads (~95%), yet unique among vertebrates [38]. These diverse morphological adaptations make frogs an ideal system in which to study modularity and integration, as they relate to ecology.

The family Myobatrachidae is an old Gondwanan lineage endemic to Australia and New Guinea with its closest relatives in South America [39]. The family currently comprises 133 described species across 21 genera, accounting for 57% of the Australian frog diversity [40]. Australia’s large landmass is characterised by a wide range of biomes and has a complex history of isolation, aridification and broad climatic changes that have had a strong impact on the evolutionary processes in its biota [41]. Myobatrachid frogs are extremely diverse in ecology (from tropical rainforest dwellers to exclusive alpine species or desert-specialists; [42]), locomotion (including excellent swimmers, jumpers, hoppers, and walkers), reproductive systems (egg deposition, calls, parental care modes, etc.; [43, 44]), and also body shape patterns [45]. Thus, they stand out as a model system to examine morphological diversification patterns on a diverse and species rich radiation across a whole continent.

We sought to address three broad questions: (i) is morphological evolution similar in different body parts, (ii) do skulls and limbs show different levels of integration?, and (iii) is morphological diversity of skulls and limbs shaped by diet, locomotion, burrowing behavior, and ecology? To do this we used 3D imaging across all genera of myobatrachids, combined with geometric morphometric analyses, to discriminate the morphological integration hypothesis and the modularity hypothesis in different structures. We used 3D data from the skull and several limb bones of the appendicular skeleton (radioulna, humerus, tibiofibula and femur), and studied their shape both as a set of distinct modules and jointly as a multi-modular integrative structure. First, we sought to quantify skull and limb shape differences across representatives of all 21 genera of myobatrachid frogs by using 3D microCT scans and geometric morphometric (GM) techniques. We then addressed three major aims. First, we tested the hypothesis that evolutionary patterns and morphological disparity are similar in the two major clades of myobatrachids across different structures. We predicted that both skull and limbs followed an evolutionary pattern typical of an adaptive radiation, and that dispersion across morphospace would be correlated with species richness, with this trend being consistent across most modules. We then determined whether there were differences in dispersion and direction of shape diversification in skulls and limbs, and whether morphological evolution acts independently in each module, or if there was some integration across different structures. We predicted a high degree of morphological integration, especially among limb modules, due to selective pressures derived from environmental correlates and associated adaptations such as burrowing behavior and locomotion. Finally, we tested for relationships between morphology and burrowing, locomotion, and environment. We predicted that form would be correlated with function, i.e. ecology, locomotion, and burrowing behavior would have been key drivers in shaping morphological evolution on the limbs, whereas head shape would be more phylogenetically conserved due to a lower functional pressure imposed by the environment.


Study samples and morphological data

This study is based on 41 ethanol-preserved specimens from 21 species of the Australo-Papuan myobatrachid frog radiation. Sampling covered all genera from this family, and with the exception of the monotypic Spicospina flammocaerulea where only one specimen was available, we used two representative specimens of the same species per genus as a previous study across all myobatrachid species showed high morphological conservatism within genera [45]. Species and voucher number details are presented in Additional file 1: Appendix S1. Since sexual dimorphism is known in some myobatrachid species (e.g. Adelotus brevis), we only sampled adult females in order to avoid morphological differences due to sexual dimorphism. We gathered data for burrowing behavior from several sources [42, 46] and classified each species into three categories based on the type of burrowing: (a) forward burrowers which use their forelimbs, (b) backwards burrowers which use their hindlimbs, and (c) non-burrowers. Locomotion information was gathered from Anstis (2013) and Cogger (2014), and locomotor mode categories were defined according to basic characteristics of their stride: (a) walkers are species that are strictly walkers or crawlers, (b) hoppers are species that can only hop, or hop and walk, and not jump (an average jumping distance that is less than five times their body length), and (c) jumpers/swimmers are species that can jump and/or swim (whose average jumping distance is greater than five times their body length and are proficient swimmers). Even though some genera display multiple states for burrowing and locomotor modes, the analyses were performed using the state present on the selected species. Data for habitat type or ecoregions was gathered taking into account each species’ distribution and the seven main ecoregions in Australia [42, 47]: (a) tropical and subtropical moist broadleaf forests, (b) temperate broadleaf and mixed forests, (c) tropical and subtropical grassland, savannas and shrublands, (d) temperate grasslands, savannas and shrublands, (e) montane grasslands and shrublands, (f) mediterranean forests, woodlands and shrubs, and (g) deserts and xeric shrublands. Dietary information [48,49,50,51,52,53,54,55,56,57,58] was gathered for all species in this study (except for the little-known species Spicospina flammocaerulea, for which we inferred diet from its close relatives and based on similarities in other life-history traits), which was classified into two categories: (a) generalists have multiple taxa represented in their diet, regardless of their size) and (b) specialists only feed on certain taxa (mostly termites and ants). Data on burrowing behavior, locomotion, diet, and ecoregions is summarized on Additional file 1: Table S1. All morphological data was gathered using three different X-ray micro-CT scanners, depending on the size of the individual frog: Skyscan 1174 (Bruker micro-CT) for small frogs, MicroXCT-400 (Xradia system) for intermediate sized frogs, and a custom-made double-helical x-ray micro CT scanner from the Australian National University for the larger specimens. The settings for each CT scanner were as follows: Skyscan 1174–40 kV source voltage, 800 μA source current, voxel size of 32.47 μm, 0.7° rotational step, 1.6 s exposure time, and 360° rotational angle scanning. The acquired images (angular projections) were reconstructed into a virtual stack of 2D cross-section slices using the NRecon (Skyscan) software interface. Xradia MicroXCT-400 - 50 kV, 360° rotational angle scanning, 2 s exposure, and voxel size of 49.13 μm. Acquired images were reconstructed in the MicroXCT and exported to a virtual stack of 2D cross-section slices (8-bit BMP format) using Avizo software system (version 8.0, Mercury Computer Systems, Inc., Germany). Custom-made double-helical x-ray micro CT - 80 kV, 100 μA, voxel size of 43 μm, using a 0.3 mm Al filter, 3.4 s exposure, and 0.143° rotational step, resulting in 2520 angular projections. This RAW data was also then reconstructed into 2D cross-section slices (NC format). Each stack of reconstructed images was then converted into 3D data, using the volume-rendering software Drishti [59].

Shape analyses

Skull and limbs’ bone shape differences were identified using geometric morphometric (GM) methods. We used rendering software Drishti [59] in order to digitise 3D landmarks of the skull and limb bones, and also sliding semilandmarks on limb bones (Additional files 2 and 3). We then averaged each dataset of morphometric data by species with geomorph [60], in order to allow analyses in a phylogenetic context and focus on morphological variation among genera and clades. We also performed GM analyses with all raw data sets before taking species means to ensure that interspecific variation was greater than intraspecific variation. Each data set was subjected to a generalised Procrustes sumperimposition fit with the package geomorph [60,61,62]. We performed a Principal Component Analysis (PCA) on the projected Procrustes coordinates into the tangent space for each set of morphological data. Each data set of GM data was analysed separately, but also joined, considering each long bone as a distinct module. To do so, we translated and rigidly rotated all landmarks and semi-landmarks from each data set using a newly developed Rigid Rotation equation, with the R package ShapeRotator [63]. This allowed us to set up all the different modules in the same position, angle and torsion and thus allow us to analyse different mobile structures as a whole (as modules would be in the same position relative to each other). We then analysed shape and size differences across all genera in each module and also in each different group of modules: (a) forelimbs (H + RU), (b) hindlimbs (F + TF) and limbs (H + RU + F + TF). In order to test our modularity and morphological integration hypotheses we also analysed morphological co-variation between: (a) skull and the four modules in the limbs (H + RU + F + TF), (b) co-variation between forelimbs (H + RU) and hindlimbs (F + TF), (c) whithin each limb, so between radioulna and humerus in forelimbs and in femur and tibiofibula in hindlimbs, and among different modules within the skull (Additional file 1: Table S2).

Statistical analyses

In order to investigate patterns of morphological evolution across the myobatrachid frog family we used a phylogeny for the group based on two mtDNA genes (ND2 and 12S) and two nDNA loci (Rag1 and Rhodopsin). This is the same phylogeny we used in a previous study of shape evolution in these frogs ([45]; toplogy available on dryad: We used the R package ape [64] to prune this tree to only include the species used in this study, and to produce an ultrametric tree with branch lengths approximating proportions of their total age. The resulting phylogeny was projected onto morphospace (previously obtained through PCA of the Procrustes coordinates) with geomorph [60] to visualise shape differences in a phylogenetic context for each of the data sets [11, 65, 66]. We also used thin-plate spline (TPS) deformation grids to visualise shape changes in the skull in the three dimension (TPS grids for x, y and x, z) using geomorph [60]. To test for the strength of phylogenetic signal in our shape data we calculated the K-statistic’s generalization for multivariate data (K mult; [67]) with geomorph [60] on the Procrustes-aligned coordinates for each GM data set. We considered a strong phylogenetic signal (K mult presenting values grater 1) as the null hypothesis which means that closely-related taxa would occupy similar regions in morphospace [68]. We tested which evolutionary model of phenotypic evolution best fits our data, for both the skull and the limbs (all four limbs bones) shape data sets, using the R packages geiger [69] and ouch [70] in the first five Principal Components. Since the results were not congruent among each PC, we decided to take a multi-variate approach using the R package mvMORPH [71], which allows complex model fitting in multivariate data. We tested the best fit for multiple models of morphological evolution in the first ten PCs of both the skull and the limbs data sets, and selected diet as a shift since it was found to be correlated to shape differences in both skulls and limbs. The models tested were: BM (Brownian Motion), BM two rates (based on diet), OU (Ornstein–Uhlenbeck), OU with two adaptive optima, EB (Early Burst), and twelve different models with a shift from two different processes at a given point in time in which some had independent rates on each time slice (Table 1).

Table 1 Summary statistics for the fit of models of phenotypic evolution in the multivariate shape datasets of Skull and limbs (all four limb bones’ analysed together): maximum likelihood estimate (ln L), sample-size corrected Akaike’s Information Criterion (AICc), and Delta AICc (∆AICc, difference between a model and the model with the lowest AICc)

We tested for evolutionary allometry by performing a regression of shape variation on size variation among different species in a phylogenetic context [72]. In order to test whether shape variation was correlated to burrowing behavior, locomotor mode, or ecoregion, we performed phylogenetic ANOVAs using the function procD.pgls() in geomorph [60] on Procrustes-aligned coordinates from each GM data set for diet, locomotor mode, burrowing behaviour, and ecology (bioregions). We also performed a phylogenetic ANOVA with all the factors, and factorial phylogentic ANOVAs with pairs of factors and their interactions. We performed a Mantel test using the R package vegan [73] to test whether there was an association between the species distribution in the skull and the limbs shape data sets, using a Spearman’s Rank correlation coefficient. Finally, we also tested for morphological disparity among the main four clades and subfamilies in the myobatrachid frog radiation in each GM data set, in relation to the number of genera per clade and the age of each clade. We used Procrustes variance (mean squared Procrustes distance of each genera from the mean shape of their clade) as a measure of morphological disparity which was calculated using geomorph [60]. Finally, we used two-block partial least squares (PLS) analysis in order to quantify shape co-variation between different structures, using geomorph [60]. We performed two-block PLS analyses between: a) skull and the overall limb shape (RU + H + TF + F), b) forelimbs (RU + H) and hindlimbs (TF + F), c) radioulna and humerus, and d) tibiofibula and femur. All two-block PLS analyses were performed on the Procrustes-aligned coordinates from each GM data set. We also assessed phylogenetic morphological integration between all these modules using the function phylo.integration() in geomorph [60].


Size and shape variation

Evolutionary allometry did not account for a significant amount of variance on skull shape: the multivariate regression of Procrustes-aligned coordinates (shape) on log-transformed centroid size (size) demonstrate that only 6.77% of the total shape variation is correlated to size variation (p = 0.23). Similarly, evolutionary allometry of limb bones was also low: only 4.29% of the total variance in total limb shape (RU + H + TF + F; p = 0.15) was correlated to size changes, 3.76% for forelimb shape (RU + H; p = 0.19), and it was slightly higher for hindlimb shape, with size correlates explaining 10.7% of the variance in shape (TF + F; p = 0.03). Given the small impact of size on shape variation we performed further analyses on the raw morphometric data sets without removing the allometric effects.

We depict skull shape variation and shape diversity across the four limb bones in Figs. 1 and 2, respectively. In the skull shape data set, the first five principal components (PC) accounted for 82.23% of the total variance (Additional file 1: Table S3), with PCskull 1 and PCskull 2 explaining 41.58% and 19.72% of morphological variation, respectively. The primary axis of variation (PCskull 1) corresponded to width and height of the cranium, and separated burrowing species (both forward burrowers and backward burrowers) and non-burrowing species (Fig. 3). The second axis of variation (PCskull 2) mainly corresponded to variation in the shape of the snout (from pointy to very rounded snouts), and clearly grouped the main two clades in different regions of the morphospace (Fig. 3). Cranium variation is also depicted in Fig. 1 through TPS grids of individuals that present the most extreme morphological variation from the consensus cranium shape.

Fig. 1

Dorsal view of skull diversity across all genera of myobatrachid frogs. The four maps display the distribution across Australian of each of the four main clades within the myobatrachids

Fig. 2

Shape diversity of limb bones in each genera of myobatrachid frogs: femur (F), tibiofibular (TF), humerus (H), and radioulna (RU). Branches on each genera have been collapsed while retaining information on the species richness of each genus. The legend depicts the three burrowing modes (forward, backward, and non-burrower) and locomotor modes (walker, hopper, and jumper/swimmer). Clades with only few species adapted to fossoriality have been indicated in the figure (Limnodynastes spp., Pseudophryne spp., and Uperoleia spp.)

Fig. 3

Phylomorphospace of PCA values on skull shape variation based on species means, using the R package geomorph. Each clade is depicted with a different shape, while burrowing behavior is indicated by different colouring. The two main different diet types (specialist and generalist) are also indicated by a schematic of each type and background colouring. Thin-plate spline (TPS) deformation grids are displayed to indicate extreme variation on skull shape among different species (names in bold), and only species from the outer limits of the morphospace are depicted

For radioulna (RU) shape variation, the first five PCs explained 78.06% of the variance (Additional file 1: Table S3), with PCRU 1 representing 57.99%, and being mostly correlated with arching on the diaphysis of the radioulna (ranging from extremely curved and constricted radioulnas in the medial part of the diaphysis to an almost straight radioulnas). PCRU 2 only added an additional 7.23% (Additional file 4: Figure S1a), and was correlated with the shape of the epiphysis. The first five PCs of the humerus (H) data set accounted for 76.79% of the overall variance (Additional file 1: Table S3), with PCH 1 representing 39.38%, and PCH 2 23.14%, mostly accounting for relative size of the deltoid tuberosity and robustness of the whole humerus (Additional file 4: Figure S1c). On the joined data set of RU and H, the first five PCs explained 82.6% of the total shape variability (Additional file 1: Table S3), with PCRU+H 1 accounting for 47.22% of the variance and PCRU+H 2 another 12.24%, and most of the morphological variability represented robustness of both humerus and radioulna, and the length of the radioulna relative to the humerus (Additional file 4: Figure S1e). On the hindlimb bones data sets, shape variation was mostly accounted within the first five PCs, with 94.82% of the total variance in tibiofibula (TF) and 97.45% in femur (F; Additional file 1: Table S3). The first axis of variation in the TF data set (PCTF 1) explained most of the morphological variation as it represented 62.01% of the overall variance (Additional file 1: Table S3) and was highly correlated with the robustness of the tibiofibula and the degree of constriction in the medial par of the diaphysis. PCTF 2 only added an additional 13.58% (Additional file 4: Figure S1b). On the F data set, PCF 1 explained 81.23% of the total morphological variance, while PCF 2 only added an additional 9.88% (Additional file 4: Figure S1d), and most of the morphological variance was correlated with the degree of arching in the medial part of the diaphysis. On the joined data set of TF and F, the first five PCs accounted for 94.72% of the variance (Additional file 1: Table S3), with PCTF+F 1 representing for 75.34% and PCTF+F 2 an additional 9.77% (Additional file 4: Figure S1f). In contrast with the TF and F data sets, the morphospace hindlimb shape axes (PCTF+F 1 and PCTF+F 2) were mostly correlated with the robustness of both the femur and tibiofibula, the amount of arching observed in the femur, the degree of constriction in the medial part of the diaphysis, and the length of the tibiofibula relative to the femur. Finally, in the overall limb bones shape data set (radioulna + humerus, + tibiofibula + femur), the first five PCs accounted for 91.19% of the morphological variation, with PClimbs 1 representing 45.46% of the variance and PClimbs 2 an additional 37.47% (Fig. 4). Most of the shape changes in the first two axes were associated with general robustness of all four bones, and were correlated with locomotor mode: walker species displayed the most negative values in both PClimbs 1 and PClimbs 2 and occupied distinct regions in the morphospace, while hoppers and jumper/swimmer species overlapped and usually displayed neutral or positive values in both axes.

Fig. 4

Phylomorphospace of PCA values on overall limb shape variation based on total shape variation of radioulna, humerus, and femur, and generated with geomorph. Different shapes depict each of the four clades, while colour represents locomotor mode: Walker, Hopper, and Jumper/Swimmer. Burrowing behavior is also indicated in the figure by a schematic of each type (burrower and non-burrower) and background colouring, and the signs ψ and * indicate whether the burrower is forward-burrowing (head and forelimbs first) or backward-burrowing (hindlimbs first), respectively. Outlines of overall body shape are displayed in species with the most extreme limb shape variation

Patterns of morphological evolution in heads and limbs

We found strong phylogenetic signal on the skull Procrustes-aligned coordinates with K mult values equivalent or greater than 1, and this was significant for the skull, femur, tibiofibular, and limbs (RU + H + TF + F; Additional file 1: Table S4). This means that more closely related species to resemble each other under a Brownian Motion process. The fitting of evolutionary models to univariate data (first five PCs) in both skull and the limb (all four limbs bones) datasets supported different models for each PC (Additional file 1: Table S5). Tests for the best fitting model of phenotypic evolution in multivariate data (first ten PCs) showed support for the same complex process in both skulls and limbs: a model of Early Burst followed by a Brownian Motion process with two different rates based on diet (Table 1). Morphological disparity of skull shape was quite similar in the two most species-rich clades, with Procrustes variance (Procvar) of 0.022 in Myobatrachinae, and Procvar = 0.030 in Limnodynastinae. In the limbs (RU + H + TF + F), morphological disparity was higher in Limnodynastinae (Procvar = 0.006) than Myobatrachinae (Procvar = 0.003). In forelimbs disparity was higher in Myobatrachinae (Procvar = 0.007) than Limnodynastinae (Procvar = 0.003). In each forelimb module separately, morphological disparity of the radioulna was higher in Myobatrachinae (Procvar = 0.007; Procvar = 0.002 in Limnodynastinae), and also in the humerus (Procvar = 0.005 in Myobatrachinae; Procvar = 0.002 in Limnodynastinae). Procrustes distances in both clades were equal in hindlimbs (TF + F; Procvar = 0.004), higher in the femurs of Limnodynastinae (Procvar = 0.004 in Myobatrachinae; Procvar = 0.006 in Limnodynastinae), and higher in the tibiofibula of Myobatrachinae (Procvar = 0.013 in Myobatrachinae; Procvar = 0.009 in Limnodynastinae). The Mantel test performed on the dissimilarity matrices extracted from the PC components of skull and limbs shape data sets was not significant (r = −0.048, p = 0.678), supporting the null hypothesis that there is no association between the species distribution in the skull and the limb morphospace.

Testing morphological integration and modularity hypotheses

The two-block partial least squares (PLS) analysis between skull and overall limb shape (RU + H + TF + F) indicates that there was slight morphological integration between head and all four limbs (r-PLS = 0.685, p = 0.011; r-PLS = 0.694, p = 0.013 after phylogenetic correction), suggesting semi-independent morphological evolution. However, this result does not hold when we look at the relationship between the head and the fore and hindlimbs separately: morphological co-variation was much higher when assessed independently on only head and forelimb (r-PLS = 0.923, p < 0.001; r-PLS = 0.909, p = 0.001 after phylogenetic correction), and even higher on head and hindlimb (r = 0.983, p = 0.002; r-PLS = 0.946, p = 0.018 after phylogenetic correction). Morphological integration between forelimbs (RU + H) and hindlimbs (TF + F) was moderate (r-PLS = 0.767, p < 0.001), but it was much higher after correcting for phylogenetic effects (r-PLS = 0.897, p = 0.001) Shape co-variation between the two modules in hindlimbs (F + TF) was extremely high (r-PLS = 0.968, p < 0.001), even after considering phylogenetic correlates (r-PLS = 0.976, p = 0.001). Similarly, the two-block PLS on the forelimbs was also high, supporting strong morphological integration between humerus and radioulna (r-PLS = 0.925, p < 0.001), even after phylogenetic correction (r-PLS = 0.932, p < 0.001). Thus, these results suggest that selective pressures acted on the two modules of hindlimbs (F + TF) and forelimbs as if it was a single integrative structure, but there was certain degree of independence between fore and hindlimbs. All the two-block PLS analyses among different modular partitions within the skull (both raw and taking phylogenetic relationships into account) displayed high levels of integration (Additional file 1: Table S5 and Additional file 5: Figure S2), pointing out that morphological features in the different substructures within the skull have evolved in concert.

Ecology, locomotion and burrowing behaviour

Phylogenetic ANOVAs performed on Procrustes-aligned coordinates of the skull data set were statistically significant for diet (F20,1 = 6.058, p = 0.001). Similarly, they were also significant for burrowing (F20,2 = 2.806, p = 0.021), and locomotor modes (F20,2 = 3.208, p = 0.001). Conversely, they were not significant for the broad eco-regions based on Australian biomes (F20,4 = 1.233, p = 0.235). In the phylogenetic ANOVA with the three significant factors (diet + burrowing + locomotion), only diet was significant (F20,1 = 3.184, p = 0.005). In the factorial phylogenetic ANOVA between burrowing and locomotion, neither the factors (F20,2 = 3.111, p = 0.157 and F20,2 = 1.952, p = 0.158, respectively) nor the interaction (F20,1 = 1.053, p = 0.316) were significant. In the factorial phylogenetic ANOVA between burrowing and diet, both diet (F20,1 = 3.368, p = 0.003) and the interaction (F20,1 = 2.373, p = 0.006) were significant, whereas burrowing was not (F20,2 = 1.551, p = 0.074). Finally, in the factorial phylogenetic ANOVA between diet and locomotion, both factors (F20,1 = 7.629, p = 0.001 and F20,2 = 2.855, p = 0.018, respectively) and the interaction (F20,1 = 2.304, p = 0.003) were significant.

On the combined limb GM data set (RU + H + TF + F) phylogenetic ANOVAs, burrowing (F20,2 = 3.113, p = 0.028), locomotion (F20,2 = 2.848, p = 0.012), and diet (F20,2 = 3.219, p = 0.006) had a significant effect on overall limb shape, whereas ecorregions (F20,4 = 1.086, p = 0.404) did not. In the phylogenetic ANOVA with combined factors of burrowing + locomotion + diet on the combined limb data set, none of the factors were significant (F20,2 = 3.265, p = 0.167; F20,2 = 1.502, p = 0.314; F20,1 = 0.877, p = 0.388; respectively). Similarly, in the factorial phylogenetic ANOVA between burrowing and locomotion, neither the factors (F20,2 = 3.161, p = 0.186 and F20,2 = 1.454, p = 0.344, respectively) nor the interaction (F20,1 = 0.371, p = 0.818) were significant. The factorial ANOVA between burrowing and diet, and the factorial phylogenetic ANOVA between diet and locomotion were also not significant. These results were slightly different when looking at forelimbs and hindlimbs data sets separately. On the forelimbs GM data set (RU + H), burrowing (F20,2 = 3.8343, p = 0.003) and diet (F20,1 = 4.383, p = 0.002) were significant, whereas locomotor mode (F20,2 = 1.310, p = 0.196) and biome were not significant (F20,4 = 0.6608, p = 0.271). In the phylogenetic ANOVA with the three factors (diet + burrowing + locomotion), only burrowing was significant (F20,2 = 5.543, p = 0.012). In the factorial ANOVA between burrowing and locomotion, only burrowing (F20,2 = 4.453, p = 0.021) was significant. In the factorial ANOVA between diet and burrowing, both factors were significant (F20,1 = 5.534, p = 0.004 and F20,2 = 3.448, p = 0.018, respectively) but the interaction was not (F20,1 = 1.095, p = 0.292). Finally, in the factorial ANOVA between diet and locomotion, only diet (F20,1 = 4.202, p = 0.012) was significant. Finally, on the hindlimbs GM data set, burrowing (F20,2 = 5.177, p = 0.013) was also significant, whereas locomotor mode (F20,2 = 1.316, p = 0.251), diet (F20,1 = 0.881, p = 0.252), and biome (F20,4 = 0.708, p = 0.687) were not. In the phylogenetic ANOVA with the three factors (burrowing + locomotion + diet) on the hindlimb GM data set, only burrowing was significant (F20,2 = 6.361, p = 0.038). In the factorial ANOVA between burrowing and locomotion, only burrowing was significant (F20,2 = 7.084, p = 0.027). None of the factors nor the interactions were significant in the factorial ANOVAs between burrowing and diet, and diet and locomotion.


We evaluated morphological differences in skulls and limb bones on representative species from all 21 genera of Australian myobatrachid frogs, using a 3D geometric morphometric approach on multiple structures. With this method we were able to focus on the tempo and mode of morphological evolution in this old Gondwanan radiation by asking three main questions: (1) whether morphological evolutionary patterns are similar for different structures, (2) if the amount and direction of morphological change differs for each structure and clade, and (3) if morphological evolution is correlated to functional traits such as locomotion, burrowing, or diet. We found that both head and limbs followed a complex evolutionary pattern typical of adaptive radiation, followed by a Brownian Motion process. Nevertheless, there was a low level of morphological integration between the skull and the limbs and there were significant differences in the mode of morphological evolution between the head and limbs. Skull morphology was phylogenetically conserved and correlated to diet, whereas limb morphology was more labile within clades and appeared to be shaped by diet, burrowing behavior and locomotion. Morphological differences among different limb modules suggest co-variation and strong morphological integration due to selection and functional constraints imposed by burrowing and locomotion. Our results illustrate how morphological diversification in animal clades can follow complex processes, entailing selective pressures from the environment as well as multiple trait covariance with varying degrees of independence across different structures. We discuss each of these topics in turn, and suggest that accurately quantifying shape diversity across multiple structures is crucial in order to understand complex evolutionary processes.

We showed that different phylogenetic clades were separated in skull morphospace, suggesting an early diversification of head shape in myobatrachid frogs, which was supported by an Early Burst model of phenotypic evolution followed by a Brownian Motion process. The majority of skull differences were correlated with fenestration: the subfamily Limnodynastinae displayed bigger and rounder orbits, and robust sphenethmoids and parasphenoids, while species from the Myobatrachinae subfamily generally showed more elongate orbits and larger antorbital fenestrae. The two other major clades, Rheobatrachus (comprising the two extinct gastric-brooding frog species) and Mixophyes spp. (8 extant species of barred frogs), displayed skull shapes that were intermediate to Limnodynastinae and Myobatrachinae. This pattern of early morphological diversification suggests that occupancy of new morphospace regions by ancestral lineages of myobatrachids could have been associated with major ecological niche filling processes that are typical of diversifying lineages [11, 74, 75]. Analogous with skull shape diversification, morphological evolution in the limbs was best explained by a complex model of an Early Burst process, followed by Brownian Motion. This was unexpected, as myobatrachid frogs are an old Gondwanan adaptive radiation, displaying an exceptionally high degree of ecological, behavioural and morphological diversity across the whole Australian continent. However, phylogenetic non-independence of highly dimensional data, such as 3D GM data, involves some potential pitfalls when inferring complex evolutionary processes, as exposed by Uyeda et al. [76]. For example, Early Burst processes could arise as an artefact from discretising highly dimensional data sets and examining a relatively small sample from multivariate patterns. Thus, caution should be used when inferring evolutionary processes on high effective dimensionality, and our initial results on myobatrachid clades would most likely benefit from a more extensive sampling within genera.

Myobatrachid frogs have experienced several major geological and climatic processes that would have affected diversification and ecological and morphological evolution [41]. Our PCA analyses of the humerus, radioulna, and whole forelimb (H + RU), distributed forward burrowers and most walkers in one broad region of the morphospace, with jumpers/swimmers, backward and non-burrowers, and one walker grouped together on the opposite side of the morphospace. The most extreme forelimb shape was exhibited by forward burrowers that displayed a stronger and more robust humerus, with larger lateral epicondyles, extremely robust radioulnas with large olecranons and a conspicuous longitudinal groove between the radius and ulna. In contrast, good jumpers or swimmer species from wet environments, such as Lechriodus sp. and the extinct Rheobatrachus spp., generally displayed slender forelimb bones with less pronounced arching, and a faint longitudinal groove in the radioulna. For the hindlimbs (F + TF), shape diversity was mostly strongly correlated with burrowing behavior (both forward burrowers and clades in which all species are backward burrowers). Both the femur and tibiofibula were shorter and thicker in burrowing species, and displayed pronounced arching and tuberosities (such as the third trochanter) to facilitate muscle attachments.

We found that diet, burrowing, and locomotion played an important role in shaping morphological diversification of myobatrachid frogs. Skull morphology was associated with diet, with ant and termite specialist feeders displaying shorter snouts than generalist species. Several other taxa, including lizards [77], crocodiles [78], mammals [22, 23], and turtles [79] also show clear associations between diet and skull shape. We also found a strong correlation between skull shape and functional traits, such as burrowing and locomotion. That is not an unexpected result, while ecotype, habitat, and other environmental and climate has been found to not have an impact in skull shape diversification in some clades [80], it can also greatly influence head shape in some clades [79, 81, 82], it can also have no impact in others.

Shape diversification of limb bones was not as strongly correlated with phylogenetic history, and instead, diet, locomotor type and burrowing behavior seemed to be important contributors to the morphological variation observed among species. Even though each limb bone displayed slight differences in their shape diversification and its correlation with different ecological variables, they did not differ substantially overall, probably due to their high morphological integration. There was, however, certain degree of independence between fore and hindlimbs, mostly due to functional requirements. Both fore and hindlimbs were correlated to burrowing behavior, but only forelimbs were associated with dietary requirements. The fact that locomotion was strongly correlated with the overall limb shape but not each module or fore and hind limbs independently is not surprising, as fore-to-hind-limb ratios have been proved to be important in explaining locomotor abilities in different frog clades [37]. While most frogs and toads explosive jumping energy is produced by the hind-limbs, they land on their adducted forelimbs, which play a critical role in locomotion by determining the landing and stabilizing actions that enable the next jumping phase [8385]. Thus, our results suggest that limb shape might have evolved as a response to locomotion constraints imposed by different structural habitats, which would constrain the locomotor modes. This concurs with results found in other amphibian clades, where variation in habitat use and locomotor behaviour seem to correlate with particular ratios between fore and hindlimb lengths [26, 27, 37]. The same trend also is noticeable in other vertebrate groups, such as phrynosomatid lizards [86], anoles [87], sauropods [88], and carnivorous mammals [89, 90], in which ecotype or locomotor type is correlated with limb morphology or distinct proportions between fore and hindlimbs.

The study of locomotion is fundamental to understanding animal biology, as it links morphology with the use of different environments through navigation, feeding, and escape from predators [91]. In addition, factors such as ecology and some less-conspicuous behavioural aspects could also contribute to morphological evolution, making inferences about evolutionary history difficult [24]. Because multiple variables can create selective pressures in different directions on phenotypic traits, their interactions could potentially lead to trade-offs. For example, morphological optimisation for burrowing creates opposing pressures from optimisation for jumping, due to discordances in functional morphological requirements for each behaviour [38]. Although myobatrachid frogs generally display phylogenetic conservatism in morphology, burrowing behavior and other ecological correlates still appear to have a strong effect on limb shape. Morphological adaptations in forward burrowers are primarily associated with forelimb bones in both amphibians and other fossorial vertebrates [92]. Similarly, despite not being found in any other vertebrate, backward burrowing represents 95% of all burrowing types in anurans [34]. The evolution of both forward and backward burrowing likely led to reduced length and increased robustness of fore and hindlimbs respectively, which would almost certainly have resulted in reduced locomotor abilities [30, 38]. This trend towards shorter and more robust limbs in burrowing anurans is likely also a beneficial adaptation to arid environments. Amphibians have adapted to a wide range of extreme climatic conditions, despite experiencing more constraints than any other terrestrial vertebrate due to rapid evaporative water loss through their permeable skin [34]. By reducing limb length, total surface area of the body can also be reduced and with it, evaporative water loss.

Despite high overall morphological disparity among different myobatrachid genera [45], some structures (e.g. limbs) displayed morphological integration and co-variation leading to convergent phenotypes, while other structures (e.g. skulls) followed semi-independent evolutionary processes. Despite the low integration between skull and all post-cranial modules, hind- and forelimbs were more tightly correlated to the skull when assessed independently, especially the hindlimbs. These results could be due to a certain degree of integration between head and postcranial modules, which could follow different directions in the morphospace for each limbs module, resulting in semi-independent pattern of morphological evolution of the head versus the rest of the body. Our results, therefore, support the modularity or semi-independent hypothesis when looking at morphological evolution between skulls and limbs, but favours the morphological integration hypothesis for shape diversification within different limb modules, or the different substructures within the skull. Thus, while high evolutionary lability experienced by limbs is a result of selective pressures from the environment, skulls instead display relatively high phylogenetic conservatism. This suggests that morphological diversification might have occurred rapidly quite early in the myobatrachid frog radiation, followed by a decrease in shape disparity, which is conspicuous through the different areas of skull morphospace. Head shape in anurans appears to have undergone extreme morphological change very early in the evolutionary history of modern amphibians, which is especially conspicuous through a substantial widening of the skull and orbits, and enlargement of fenestrae [76, 93]. Moreover, strong phylogenetic structure on skull shape is not unusual among other amphibian groups older than 50 MY (e.g. caecilians [75]), in contrast to younger vertebrate radiations that typically display greater morphological disparity, with weaker phylogenetic signal [94].

Phylogenetic conservatism and morphological diversification in functional traits can provide insight into evolutionary processes [24], but the interplay between different potential drivers of adaptation can blur the link between form and function. For example, limb morphology might appear strongly correlated with locomotion type, but habitat use or burrowing behaviour might be equally important correlates. In this way adaptive traits often cannot (and should not) be explained by just one adaptive process. Morphological integration or modularity also can affect the accuracy of evolutionary inferences on adaptation to certain ecological, locomotor or behavioural factors [95]. Furthermore, rates of phenotypic evolution can be correlated with species diversification rates within clades, as morphological traits typically have slower evolutionary rates than other traits such as behaviour [67, 96]. Moreover, closely related clades might display unequal magnitudes of morphological change, thus hindering or boosting apparent morphological diversification, especially early in their evolutionary history [11].


Our study is the first to accurately identify evolutionary processes that drive morphological diversity in the context of modularity and morphological integration of several structures in an old adaptive radiation and across a whole continent. Our results highlight how form is usually tightly linked to function, and that different structures can evolve semi-independently, while in other modules morphological evolution is tightly coupled. There was strong morphological co-variation among different modules in the limbs due to strong selective pressures from the environment and functional trade-offs (e.g. burrowing and locomotion), whereas skull shape was correlated to diet, and yield a pattern of very early morphological diversification followed by strong phylogenetic conservatism. Our results also show that even when different structures evolve following the same evolutionary models, patterns of morphological diversification can be drastically different. The complex interplay between selective pressures and different levels of co-varying morphological evolution makes it harder to accurately identify processes that drive clade diversification and infer their evolutionary history. Thus, we highlight the importance of accurately assessing morphological evolution in multiple structures in order to properly understand complex evolutionary processes that generate the phenotypic diversity we see today.


  1. 1.

    Russell ES. Form and function. A contribution to the history of animal morphology. 1917;

  2. 2.

    La Barbera M. Analysing body size as a factor in ecology and evolution. Annu Rev Ecol Syst. 1989;20:97–117.

    Article  Google Scholar 

  3. 3.

    Blackburn TM, Gaston KJ. Animal body size distributions: patterns, mechanisms and implications. Trends Ecol Evol. 1994;9:471–4.

    CAS  Article  PubMed  Google Scholar 

  4. 4.

    Collar DC, Near TJ, Wainwright PC. Comparative analysis of morphological diversity: does disparity accumulate at the same rate in two lineages of centrarchid fishes? Evolution. 2005;59:1783–94.

    Article  PubMed  Google Scholar 

  5. 5.

    Ricklefs RE. Community diversity: relative roles of local and regional processes. Science. 1987;235:167–71.

    CAS  Article  PubMed  Google Scholar 

  6. 6.

    Wainwright PC. Ecomorphology: experimental functional anatomy for ecological problems. Am Zool. 1991;31:680–93.

    Article  Google Scholar 

  7. 7.

    Wainwright PC, Reilly SM. Ecological morphology: integrative organismal biology [internet]: University of Chicago Press; 1994.

  8. 8.

    Losos JB. The evolution of form and function: morphology and locomotor performance in West Indian Anolis lizards. Evolution. 1990;44:1189–203.

    Article  PubMed  Google Scholar 

  9. 9.

    Foote M. The evolution of morphological diversity. Annu. Rev. Ecol. Syst. 1997;28:129–52.

    Google Scholar 

  10. 10.

    Crisp MD, Arroyo MTK, Cook LG, Gandolfo MA, Jordan GJ, McGlone MS, et al. Phylogenetic biome conservatism on a global scale. Nature. 2009;458:754–6.

    CAS  Article  PubMed  Google Scholar 

  11. 11.

    Sidlauskas B. Continuous and arrested morphological diversification in sister clades of characiform fishes: a phylomorphospace approach. Evolution. 2008;62:3135–56.

    Article  PubMed  Google Scholar 

  12. 12.

    Thacker CE. Species and shape diversification are inversely correlated among gobies and cardinalfishes (Teleostei: Gobiiformes). Org Divers Evol. 2014;14:419–36.

    Article  Google Scholar 

  13. 13.

    Claramunt S. Discovering exceptional diversifications at continental scales: the case of the endemic families of neotropical suboscine passerines. Evolution. 2010;64:2004–19.

    PubMed  Google Scholar 

  14. 14.

    Brusatte SL, Nesbitt SJ, Irmis RB, Butler RJ, Benton MJ, Norell MA. The origin and early radiation of dinosaurs. Earth-Science Rev. 2010;101:68–100.

    Article  Google Scholar 

  15. 15.

    Klingenberg CP. Morphological integration and developmental modularity. Annu Rev Ecol Evol Syst. 2008;39:115–32.

    Article  Google Scholar 

  16. 16.

    Mitteroecker P, Bookstein F. The conceptual and statistical relationship between modularity and morphological integration. Syst Biol. 2007;56:818–36.

    Article  PubMed  Google Scholar 

  17. 17.

    Gatesy SM, Dial KP. Locomotor modules and the evolution of flight. Evolution. 1996;50:331–40.

    Article  PubMed  Google Scholar 

  18. 18.

    Feilich KL. Correlated evolution of body and fin morphology in the cichlid fishes. Evolution. 2016;70:2247–67.

    Article  PubMed  Google Scholar 

  19. 19.

    Hansen TF. Is modularity necessary for evolvability? Biosystems. 2003;69:83–94.

    Article  PubMed  Google Scholar 

  20. 20.

    Raff RA. The shape of life: genes, development, and the evolution of animal form: University of Chicago Press; 2012.

  21. 21.

    Fruciano C, Franchini P, Meyer A. Resampling-based approaches to study variation in morphological modularity. PLoS One. 2013;8:e69376.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  22. 22.

    Goswami A, Polly PD. The influence of modularity on cranial morphological disparity in Carnivora and primates (Mammalia). PLoS One. 2010;5:e9517.

    Article  PubMed  PubMed Central  Google Scholar 

  23. 23.

    Goswami A. Notes and comments cranial modularity shifts during mammalian. Evolution. 2012;168:270–80.

    Google Scholar 

  24. 24.

    Ricklefs RE, Miles DB. Ecological and evolutionary inferences from morphology: an ecological perspective. in Ecological morphology: integrative organismal biology. University of Chicago Press; 1994. p. 13–41.

  25. 25.

    Zaaf A, Van Damme R. Limb proportions in climbing and ground-dwelling geckos (Lepidosauria, Gekkonidae): a phylogenetically informed analysis. Zoomorphology. 2001;121:45–53.

    Article  Google Scholar 

  26. 26.

    Moen DS, Irschick DJ, Wiens JJ. Evolutionary conservatism and convergence both lead to striking similarity in ecology, morphology and performance across continents in frogs. Proc Biol Sci. 2013;280:2013–156.

    Article  Google Scholar 

  27. 27.

    Enriquez-Urzelai U, Montori A, Llorente GA, Kaliontzopoulou A. Locomotor mode and the evolution of the Hindlimb in western Mediterranean anurans. Evol Biol. 2015;42:199–209.

    Article  Google Scholar 

  28. 28.

    Cornette R, Baylac M, Souter T, Herrel A. Does shape co-variation between the skull and the mandible have functional consequences? A 3D approach for a 3D problem. J Anat. 2013;223:329–36.

    Article  PubMed  PubMed Central  Google Scholar 

  29. 29.

    Wagner GP, Pavlicev M, Cheverud JM. The road to modularity. Nat Rev Genet. 2007;8:921–31.

    CAS  Article  PubMed  Google Scholar 

  30. 30.

    Hall BK. Fins into limbs: evolution, development, and transformation: University of Chicago Press; 2008.

  31. 31.

    Shubin NH, Jenkins FA. An early Jurassic jumping frog. Nature. 1995;377:49–52.

    CAS  Article  Google Scholar 

  32. 32.

    Jenkins F A., Shubin NH. Prosalirus bitis and the anuran caudopelvic mechanism. J Vertebr Paleontol 1998;18:495–510.

  33. 33.

    Reilly S, Essner R, Wren S, Easton L, Bishop PJ. Movement patterns in leiopelmatid frogs: insights into the locomotor repertoire of basal anurans. Behav Process. 2015;121:43–53.

    Article  Google Scholar 

  34. 34.

    Wells KD. The ecology and behavior of amphibians: University of Chicago Press; 2010.

  35. 35.

    Emerson SB, Koehl MAR. The interaction of behavioral and morphological change in the evolution of a novel locomotor type: “flying” frogs. Evolution. 1990;44:1931–46.

    PubMed  Google Scholar 

  36. 36.

    Wilson RS, James RS, Van Damme R. Trade-offs between speed and endurance in the frog Xenopus laevis: a multi-level approach. J Exp Biol. 2002;205:1145–52.

    PubMed  Google Scholar 

  37. 37.

    Vidal-García M, Keogh JS. Convergent evolution across the Australian continent: ecotype diversification drives morphological convergence in two distantly related clades of Australian frogs. J Evol Biol. 2015;28:2136–51.

    Article  PubMed  Google Scholar 

  38. 38.

    Emerson SB. Burrowing in frogs. J Morphol. 1976;149:437–58.

    Article  Google Scholar 

  39. 39.

    Pyron RA. Biogeographic analysis reveals ancient continental Vicariance and recent oceanic dispersal in amphibians. Syst Biol. 2014;63:779–97.

    Article  PubMed  Google Scholar 

  40. 40.

    Frost DR. Amphibian Species of the World: an Online Reference. Version 6.0. [cited. Oct 15. 2016; Available from:

  41. 41.

    Byrne M, Yeates DK, Joseph L, Kearney M, Bowler J, Williams MAJ, et al. Birth of a biome: insights into the assembly and maintenance of the Australian arid zone biota. Mol Ecol. 2008;17:4398–417.

    CAS  Article  PubMed  Google Scholar 

  42. 42.

    Anstis M. Tadpoles and frogs of Australia. Sydney, NSW: New Holland Publishing Pty Ltd; 2013.

    Google Scholar 

  43. 43.

    Roberts J, Standish R, Byrne P, Doughty P. Synchronous polyandry and multiple paternity in the frog Crinia georgiana (Anura: Myobatrachidae). Anim Behav. 1999;57:721–6.

    CAS  Article  PubMed  Google Scholar 

  44. 44.

    Byrne PG, Roberts JD, Simmons LW. Sperm competition selects for increased testes mass in Australian frogs. J Evol Biol. 2002;15:347–55.

    Article  Google Scholar 

  45. 45.

    Vidal-García M, Byrne PG, Roberts JD, Keogh JS. The role of phylogeny and ecology in shaping morphology in 21 genera and 127 species of Australo-Papuan myobatrachid frogs. J Evol Biol. 2014;27:181–92.

    Article  PubMed  Google Scholar 

  46. 46.

    Cogger H. Reptiles and amphibians of Australia. 7th ed. Collingwood, VIC: Csiro Publishing; 2014.

    Google Scholar 

  47. 47.

    National Reserve System. National Reserve System (NRS). 2016. Available from:

  48. 48.

    Harrison L. Notes on some western Australian frogs, with descriptions of new species. Rec Aust Museum. 1927;15:277–87.

    Article  Google Scholar 

  49. 49.

    Calaby JH. The food habits of the frog, Myobatrachus gouldii (Gray). West Aust Nat. 1956;5:93–4.

    Google Scholar 

  50. 50.

    Calaby JH. A note on the food of Australian desert frogs. West. Aust. Nat. 1960;7:79–80.

    Google Scholar 

  51. 51.

    Lee AK. Studies in Australian amphibia II..Taxonomy, ecology and evolution of the genus Heleioporus Gray (Anura : Leptodactylidae). Aust. J. Zool 1967;15:367–439.

  52. 52.

    Tyler MJ, Roberts JD, Davies M. Field observations on Arenophryne rotunda Tyler, a Leptodactylid frog of coastal Sandhills. Aust Wildl Res. 1980;7:295–304.

    Article  Google Scholar 

  53. 53.

    Mac Nally RC. Trophic relationships of two sympatric species of Ranidella (Anura). Herpetologica. 1983:130–40.

  54. 54.

    Winter J, McDonald R. Eungella, the land of cloud. Aust Nat Hist. 1986;22:39–43.

    Google Scholar 

  55. 55.

    Cappo M. Frogs as predators of organisms of aquatic origin in the Magela Creek system. Northern Territory. Thesis: University of Adelaide, Department of Zoology; 1987.

    Google Scholar 

  56. 56.

    Katsikaros K, Shine R. Sexual dimorphism in the tusked frog, Adelotus brevis (Anura:Myobatrachidae): the roles of natural and sexual selection. Biol J Linn Soc. 1997;60:39–51.

    Google Scholar 

  57. 57.

    Lima AP, Magnusson WE, Williams DG. Differences in diet among frogs and lizards coexisting in subtropical forests of Australia. J Herpetol. 2000;34:40–6.

    Article  Google Scholar 

  58. 58.

    Mahony M, Clulow S, Roberts JD. Personal communication. 2016.

    Google Scholar 

  59. 59.

    Limaye A. Drishti: a volume exploration and presentation tool. SPIE 8506, Dev. X-Ray Tomogr. 2012;8506:85060X.

  60. 60.

    Adams DC, Otárola-Castillo E. geomorph : an R package for the collection and analysis of geometric morphometric shape data. 2013;63:685–697.

  61. 61.

    Klingenberg CP, Gidaszewski N. Testing and quantifying phylogenetic signals and homoplasy in morphometric data. Syst Biol. 2010;59:245–61.

    CAS  Article  PubMed  Google Scholar 

  62. 62.

    Rohlf F, Slice D. Extensions of the Procrustes method for the optimal superimposition of landmarks. Syst Biol. 1990;39:40–59.

    Google Scholar 

  63. 63.

    Vidal-Garcia M, Bandara L, Keogh JS. ShapeRotator: an R package for standardised rigid rotations of articulated Three-Dimensional structures with application for geometric morphometrics. bioRxiv. 2017.

  64. 64.

    Paradis E, Claude J, Strimmer K. APE: analyses of Phylogenetics and evolution in R language. Bioinformatics. 2004;20:289–90.

    CAS  Article  PubMed  Google Scholar 

  65. 65.

    Maddison WP. Squared-change parsimony reconstructions of ancestral states for continuous-valued characters on a phylogenetic tree. 2013;40:304–314.

  66. 66.

    Rohlf FJ. Geometric morphometrics and phylogeny. In: MacLeod N, Forey PL, editors. Morphol. shape phylogeny. 2002. p. 175–193.

  67. 67.

    Adams DC. A generalized K statistic for estimating phylogenetic signal from shape and other high-dimensional multivariate data. Syst Biol. 2014;63:685–97.

    Article  PubMed  Google Scholar 

  68. 68.

    Blomberg SP, Garland T, Ives AR. Testing for phylogenetic signal in comparative data: behavioral traits are more labile. Evolution. 2003;57:717–45.

    Article  PubMed  Google Scholar 

  69. 69.

    Harmon LJ, Weir JT, Brock CD, Glor RE, Challenger W. GEIGER: investigating evolutionary radiations. Bioinformatics. 2007;24(1):129–31.

  70. 70.

    King AA, Butler MA. ouch: Ornstein-Uhlenbeck models for phylogenetic comparative hypotheses (R package). 2009. Available from:

  71. 71.

    Clavel J, Escarguel G, Merceron G. mvMORPH: an R package for fitting multivariate evolutionary models to morphometric data. Methods Ecol. Evolution. 2015;6:1311–9.

    Google Scholar 

  72. 72.

    Klingenberg CP, Ekau W. A combined morphometric and phylogenetic analysis of an ecomorphological trend: pelagization in Antarctic fishes (Perciformes: Nototheniidae). Biol J Linn Soc. 1996;59:143–77.

    Article  Google Scholar 

  73. 73.

    Dixon P. Vegan, a package of R functions for community ecology. J Veg Sci. 2003;14:927–30.

    Article  Google Scholar 

  74. 74.

    Schluter D. Ecological character displacement in adaptive radiation. Am Nat. 2000;156:S4–S16.

    Article  Google Scholar 

  75. 75.

    Sherratt E, Gower DJ, Peter C, Mark K. Evolution of cranial shape in caecilians (Amphibia : Gymnophiona). Evol Biol. 2014;4:528–45.

    Article  Google Scholar 

  76. 76.

    Uyeda JC, Caetano DS, Pennell MW. Comparative analysis of principal components can be misleading. Syst Biol. 2015;64(4):677–89.

  77. 77.

    Stayton CT. Morphological Evolution of the Lizard Skull : A Geometric Morphometrics Survey. 2005;59:47–59.

    Google Scholar 

  78. 78.

    Pierce SE, Angielczyk KD, Rayfield EJ. Morphospace occupation in thalattosuchian crocodylomorphs : skull shape variation, species delineation and temporal patterns. Palaeontology. 2009;52:1057–97.

    Article  Google Scholar 

  79. 79.

    Claude J, Pritchard P, Tong H, Paradis E, Auffray J-C. Ecological correlates and evolutionary divergence in the skull of turtles: a geometric morphometric assessment. Syst Biol. 2004;53:933–48.

    Article  PubMed  Google Scholar 

  80. 80.

    Kohlsdorf T, Grizante MB, Navas CA, Herrel A. Head shape evolution in Tropidurinae lizards: does locomotion constrain diet? J Evol Biol. 2008;21:781–90.

    CAS  Article  PubMed  Google Scholar 

  81. 81.

    Kaliontzopoulou A, Carretero MA, Llorente GA. Intraspecific ecomorphological variation : linear and geometric morphometrics reveal habitat-related patterns within Podarcis bocagei wall lizards. J Evol Biol. 2010;23:1234–44.

    Article  PubMed  Google Scholar 

  82. 82.

    Kaliontzopoulou A, Adams DC, van der Meijden A, Perera A. Carretero M a. Relationships between head morphology, bite performance and ecology in two species of Podarcis wall lizards. Evol. Ecol. 2012;26:825–45.

    Google Scholar 

  83. 83.

    Nauwelaerts S, Stamhuis E, Aerts P. Swimming and jumping in a semi-aquatic frog. Anim Biol. 2005;55:3–15.

    Article  Google Scholar 

  84. 84.

    Essner RL, Suffian DJ, Bishop PJ, Reilly SM. Landing in basal frogs: evidence of saltational patterns in the evolution of anuran locomotion. Naturwissenschaften. 2010;97:935–9.

    CAS  Article  PubMed  Google Scholar 

  85. 85.

    Gillis GB, Akella T, Gunaratne R. Do toads have a jump on how far they hop? Pre-landing activity timing and intensity in forelimb muscles of hopping Bufo marinus. Biol Lett. 2010;6:486–9.

    Article  PubMed  PubMed Central  Google Scholar 

  86. 86.

    Herrel A, Meyers JJ, Vanhooydonck B. Relations between microhabitat use and limb shape in phrynosomatid lizards. Biol J Linn Soc. 2002;77:149–63.

    Article  Google Scholar 

  87. 87.

    Losos JB. Ecomorphology, performance capability, and scaling of West Indian Anolis lizards: an evolutionary analysis. Ecol Monogr. 1990;60:369–88.

    Article  Google Scholar 

  88. 88.

    Bonnan MF. Morphometric analysis of humerus and femur shape in Morrison sauropods: implications for functional morphology and paleobiology. Paleobiology. 2004;30:444–70.

    Article  Google Scholar 

  89. 89.

    Ercoli MD, Prevosti FJ, Álvarez A. Form and function within a phylogenetic framework : locomotory habits of extant predators and some Miocene Sparassodonta (Metatheria). Zool J Linnean Soc. 2012;165:224–51.

    Article  Google Scholar 

  90. 90.

    Fabre A-C, Cornette R, Goswami A, Peigné S. Do constraints associated with the locomotor habitat drive the evolution of forelimb shape? A case study in musteloid carnivorans. J Anat. 2015;226:596–610.

    Article  PubMed  PubMed Central  Google Scholar 

  91. 91.

    Gray J. How animals move. Harmondsworth: Penguin Books; 1959.

    Google Scholar 

  92. 92.

    Piras P, Sansalone G, Teresi L, Kotsakis T, Colangelo P. Loy a. Testing convergent and parallel adaptations in talpids humeral mechanical performance by means of geometric morphometrics and finite element analysis. J. Morphology. 2012;273:696–711.

    CAS  Article  Google Scholar 

  93. 93.

    Sigurdsen T, Bolt JR. The lower Permian amphibamid Doleserpeton (Temnospondyli: Dissorophoidea), the interrelationships of amphibamids, and the origin of modern amphibians. J Vertebr Paleontol. 2010;30:1360–77.

    Article  Google Scholar 

  94. 94.

    Losos JB. Adaptive radiation, ecological opportunity, and evolutionary determinism. Am Nat. 2010;175:623–39.

    Article  PubMed  Google Scholar 

  95. 95.

    Young N. Modularity and integration in the hominoid scapula. J Exp Zool B Mol Dev Evol. 2004;302:226–40.

    Article  PubMed  Google Scholar 

  96. 96.

    Harmon LJ, Losos JB, Davies TJ, Gillespie RG, Gittleman JL. Jennings WB, et al. Early bursts of body size and shape evolution are rare in comparative data. 2010:2385–96.

Download references


We are grateful to all museum curators for specimens’ loans: R. Sadlier, M. Hutchinson, P. Doughty. We thank L. Joseph, R. Palmer and M. Cawsey (ANWC); M. Foley (ACMM); and M. Turner and T. Senden (ANU) for facilitating access to CT scan facilities. We also thank M. Mahony, S. Clulow and D. Roberts for their personal communications on frog diet. We are indebted to A. Pyron, M. Pepper, I. Brennan, T. Semple, and an anonymous reviewer for useful comments on this manuscript. Thank you to E. Walsh for her drawings for the burrowing legends in Fig. 4.


All the data was gathered as part of MVG’s PhD at the Australian National University. MVG was supported by the Peter Rankin Trust for Herpetology and the SSAR. JSK was funded by the ARC.

Availability of data and materials

All data will be available in the Dryad repository upon acceptance.

Authors’ contributions

MVG and SK conceived the study. MVG collected, processed, and analysed the data, and drafted the initial version of the manuscript. Both authors read, edited, and approved the final version of the manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not applicable.

Ethics approval and consent to participate

Not applicable.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Author information



Corresponding author

Correspondence to Marta Vidal-García.

Additional files

Additional file 1:

Table S1. Summary of several ecological and behavioural traits of the myobatrachid frogs studied here, used in posterior analyses: burrowing behaviour, locomotor mode, habitat type or ecoregion, and diet type. Table S2. Summary of the different landmarks used for each module (m1, m2 or m3) in all five models of modular partitions (bimodular and trimodular) within the skull that correspond to the models displayed on Additional file 5: Figure S2. Table S3. Principal Component Analyses of shape variation for different sets of Procrustes-aligned species means, using geomorph. Table S4. Summary of phylogenetic signal tests, using geomorph (Adams & Otarola-Castillo, 2013). K 95% confidence interval for values expected under a Brownian Motion model of trait evolution = [0.799, 1.318]. Table S5. Summary statistics for the fit of models of phenotypic evolution in the first five principal components of the Skull shape dataset and the limbs shape dataset (all four limb bones together): maximum likelihood estimate (ln L), sample-size corrected Akaike’s Information Criterion (AICc), and Delta AICc (ΔAICc, difference between a model and the model with the lowest AICc). We tested the fit of the following evolutionary models: BM = Brownian Motion, EB = Early Burst, white = nonphylogenetic, OU = Ornstein-Uhlenbeck, OU2_diet = Ornstein-Uhlenbeck with two optima based on diet, OU3_loc = OU with three optima based on locomotion, and OU3_burr = OU with three optima based on burrowing behaviour. Analyses were performed in R using geiger [68] and ouch [69]. Table S6. Results from the integration.test function in geomorph (Adams & Otrola-Castillo, 2013) in order to quantify the degree of modularity between the two or three modules in each modular configuration (a-e), using the landmark coordinate data. Appendix S1. Species and specimen codes for all the individuals used in this study, by museums. (PDF 171 kb)

Additional file 2:

Video displaying the 42 landmarks used in the GM analyses of the skull. (PDF 25415 kb)

Additional file 3:

Video displaying both the landmarks and semi-landmarks used in the GM analyses of the four limb bones: radioulna, humerus, tibiofibular and femur. (PDF 6874 kb)

Additional file 4: Figure S1.

(a) Phylomorphospace of PCA values on shape variation of radioulna (RU); (b) Phylomorphospace of PCA values on shape variation of tibiofibula (TF); (c) Phylomorphospace of PCA values on shape variation of humerus (H); (d) Phylomorphospace of PCA values on shape variation of femur (F); (e) Phylomorphospace of PCA values on fore-limb shape variation (RU + H); (f) Phylomorphospace of PCA values on hind-limb shape variation (TF + F). (PDF 332 kb)

Additional file 5: Figure S2.

Modular configurations modeled for the skull with two or three different partitions, based on different evolutionary hypothesis based on biological relevant regions. The different colours depict different modules. (a) The first module (green) includes the tip of the snout and the olfactory area (premaxilla, maxilla, and nasal), as it captures a lot of morphological variation among frog species, whereas the second module (blue) includes the rest of the skull; (b) this configuration captures skull depth – the first module includes the dorsal region of the skull, and the second module captures morphological information from the ventral region; (c) this tripartite model splits the skull in three modules: snout (green), squamosal (orange, which is part of the suspensory apparatus), and the rest of the skull (blue); (d) The first module depicts the snout (green), the second includes the medial part of the skull (orange), and the third module includes the most posterior region of the skull (blue); (e) this tripartite configuration includes a first module (green) with the snout morphology, a second module (orange) that encompasses the brain region (from the sphenethemoid to the exoccipital and foramen magnum, including the frontoparietal), and a third module (blue) for the rest of the skull. (PDF 275 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Vidal-García, M., Scott Keogh, J. Phylogenetic conservatism in skulls and evolutionary lability in limbs – morphological evolution across an ancient frog radiation is shaped by diet, locomotion and burrowing. BMC Evol Biol 17, 165 (2017).

Download citation


  • Morphology
  • Modularity
  • Morphological integration
  • 3D morphology
  • Geometric morphometrics
  • Phylomorphospace