Caecilians are elongate, limbless and annulated amphibians that, with the exception of one aquatic family, all have an at least partly fossorial lifestyle. It has been suggested that caecilian evolution resulted in sturdy and compact skulls with fused bones and tight sutures, as an adaptation to their head-first burrowing habits. However, although their cranial osteology is well described, relationships between form and function remain poorly understood. In the present study, we explored the relationship between cranial shape and in vivo burrowing forces. Using micro-computed tomography (µCT) data, we performed 3D geometric morphometrics to explore whether cranial and mandibular shapes reflected patterns that might be associated with maximal push forces. The results highlight important differences in maximal push forces, with the aquatic Typhlonectes producing a lower force for a given size compared with other species. Despite substantial differences in head morphology across species, no relationship between overall skull shape and push force could be detected. Although a strong phylogenetic signal may partly obscure the results, our conclusions confirm previous studies using biomechanical models and suggest that differences in the degree of fossoriality do not appear to be driving the evolution of head shape.
Caecilians (Gymnophiona) are a monophyletic group of elongate, limbless and annulated amphibians (Channing, 2001; Pough et al., 1998; Taylor, 1968). As most species are fossorial to some degree, their ecology remains generally poorly known (Nussbaum and Wilkinson, 1989; O'Reilly, 2000; Summers and O'Reilly, 1997; Taylor, 1968). Among the 215 extant species that are classified into 10 families (AmphibiaWeb: https://amphibiaweb.org, accessed 23 July 2021; Kamei et al., 2012; Wilkinson et al., 2011), only Typhlonectidae includes aquatic and semi-aquatic species, the nine other families being terrestrial and more or less fossorial (Taylor, 1968). The cranial system plays vital roles in many activities, including feeding and respiration, in addition to housing and protecting the brain and major sensory organs (Wake, 1993). In limbless fossorial vertebrates, head-first burrowing imposes additional constraints on the cranial system (Maddin et al., 2012; O'Reilly, 2000; Wake, 1993). Indeed, in direct association with their burrowing habits, it is thought that caecilian skull evolution resulted in compact and robust skulls with bones being fused or connected through tight sutures (Nussbaum and Pfrender, 1998; Taylor, 1968; Wake, 1993; Wake and Hanken, 1982).
Despite a relatively conserved morphology and homoplasy among some clades, caecilians radiated into a variety of ecological niches (Sherratt et al., 2014; Taylor, 1968). Associated with this ecological variation, the position of their mouth varies from terminal to sub-terminal. Because the head is used for the initial substrate penetration, sub-terminal mouths are expected to be found in more active burrowers (Gans, 1973, 1974; Hohl et al., 2014; O'Reilly, 2000; Wake, 1993). Although all caecilians have reduced eyes, the orbits are completely covered with bone in some species (Mohun and Wilkinson, 2015; Wake, 1985). Two different skull morphologies are found in caecilians: zygokrotaphic and stegokrotaphic skulls (Kleinteich et al., 2012; Sherratt et al., 2014; Wake, 2003), with zygokrotaphic skulls showing an incomplete coverage of the temporal region by the squamosal, as observed in Scolecomorphidae, Typhlonectidae and Rhinatrematidae (Kleinteich et al., 2012; Nussbaum, 1977, 1985; Sherratt et al., 2014; Wilkinson and Nussbaum, 1997). In animals with stegokrotaphic skulls, the squamosal completely covers the temporal region, as observed in Siphonopidae, Indotyphlidae, Caecilidae, Chikilidae, Herpelidae and Ichthyophiidae (Sherratt et al., 2014; Wake, 2003). The Dermophiidae comprise species with zygokrotaphic and stegokrotaphic skulls (Taylor, 1968; Wake, 2003). However, although stegokrotaphic skulls were thought to be structurally better suited to resist forces encountered during head-first burrowing (Nussbaum, 1977, 1983), previous studies found no strong evidence for this hypothesis (Ducey et al., 1993; Herrel and Measey, 2010; Kleinteich et al., 2012). Moreover, some species with stegokrotaphic skulls are not excellent burrowers (Ducey et al., 1993; Herrel and Measey, 2010; Nussbaum and Pfrender, 1998; Wollenberg and Measey, 2009) and are found mostly in the leaf litter (Burger et al., 2004; Gower et al., 2004; Kupfer et al., 2005).
Burrowing is a complex behaviour that remains rather poorly understood (but see Gans, 1973; Gaymer, 1971; Hohl et al., 2014; Summers and O'Reilly, 1997; Teodecki et al., 1998). Thanks to a partial independence between the body wall and vertebral column, most caecilians are capable of moving through narrow tunnels using a combination of hydrostatic and internal concertina locomotion (Gaymer, 1971; O'Reilly et al., 1997; Summers and O'Reilly, 1997; Herrel and Measey, 2010). In addition to being able to move through existing tunnels, the ability to actively create tunnels is probably a key trait in fossorial vertebrates, as it allows them to explore their environment and find food resources and sexual partners. Although burrowing speed is likely to be a key performance trait when escaping predators, an increase in the maximal forces that can be exerted may be critical to explore a broader variety of soils.
Here, we present data on the maximal push forces exerted by caecilians from different families to test whether: (1) aquatic caecilians have a lower performance than terrestrial species, as predicted given that these animals have lost the ability to perform internal concertina locomotion (Summers and O'Reilly, 1997); and (2) phylogenetically basal rhinatrematids can produce forces equal to those observed in other families. Despite the fact that these animals are known to show skin–vertebral independence they are an early-diverging lineage and have an unspecialized head shape with a terminal mouth (Naylor and Nussbaum, 1980). We then quantify the cranial and mandibular shape using 3D geometric morphometrics and investigate the relationships between shape and in vivo burrowing forces. In accordance with studies on other limbless reptilian burrowers (i.e. skinks: Vanhooydonck et al., 2011; amphisbaenians: Hohl et al., 2017), we also predict that for caecilians, narrower, more fusiform heads would be associated with greater forces, as these shapes probably facilitate soil penetration (e.g. Measey and Herrel, 2006; Herrel and Measey, 2010; Barros et al., 2011). Alternatively, species with blunt heads may need to be able to produce greater forces to penetrate the substrate, putting a premium on the ability to generate push force.
MATERIALS AND METHODS
We quantified the shape of the cranium and mandible of 81 individuals from 28 species belonging to nine out of the 10 currently recognized families (Table 1), thus capturing a broad diversity in cranial osteology, phylogeny and ecology. Our sample was restricted to adults and included both males and females. Although a sexual dimorphism is present in caecilians (Kupfer, 2009; Maerker et al., 2016), interspecific variation largely exceeds the sex-specific variation (Sherratt et al., 2014). Specimens were primarily obtained from our personal collections and completed with specimens from museum collections (Table S1).
Micro-computed tomography (μCT) imaging
For this study, a large amount of CT scans of different species was used, as listed in Table S1. About half of these scans were performed at the Centre for X-Ray Tomography at Ghent University, Belgium (UGCT, www.ugct.ugent.be) using the HECTOR µCT scanner (Masschaele, 2013). The scanner settings were sample dependent. The tube voltage varied between 100 and 120 kV and the number of X-ray projections taken over 360 deg was typically about 2000 per scan. Additional µCT scans were obtained from the online repository Morphosource (morphosource.org), the Zoological Museum Hamburg (see Kleinteich et al., 2008a,b for scanner settings), the royal Museum of Central Africa (75 kV, 1440 projections) and the personal collection of Mark Wilkinson (100 kV, 2000 projections; see Table S1). The isotropic voxel size of all scans is listed in Table S1. All the µCT-scans were processed using both automatic thresholding and manual segmentation to reconstruct the cranium and mandible in 3D using Amira 2019.3 (Visage Imaging, San Diego, CA, USA). Using Geomagic Wrap (3D systems), surfaces were prepared by removing highly creased edges and spikes that may interfere with the placement of landmarks. Next, they were decimated to a maximum of approximately 700,000 faces to reduce computational demands without compromising details. For the mandible, only the left hemi-mandible was used. The ‘mirror’ function in Geomagic was used for the specimens where the left hemi-mandible was damaged.
In vivo burrowing forces were measured in the field for 120 specimens belonging to 13 species (Table 1). Specimens were maintained for a maximum of 24 h in large containers filled with substrate collected in the field at sites where the animals were found. After measurements, animals were released at the exact locations where they were found. Push forces were measured using a custom piezoelectric force platform (Fig. 1; Kistler Squirrel force plate, ±0.1 N, Kistler, Zurich, Switzerland) as described in Vanhooydonck et al. (2011). The force plate was mounted on a purpose-built metal base (Fig. 1A) and connected to a Kistler charge amplifier (type 9865). A Perspex block with 1 cm deep holes of different diameters was mounted on the force plate (Fig. 1B), level with the front edge. One hole with a diameter corresponding to the body diameter of the animal was filled with substrate from the container of the animal being tested. A tunnel with a diameter equal to the maximum body width of the animal was positioned on the metal base in front of (but not touching) the soil-filled hole in the Perspex block. Then, a caecilian was introduced into the tunnel and allowed to move towards the Perspex block. Next, the animal was stimulated to push into the soil-filled hole by gently tapping the end of its tail region. Forces were recorded in three dimensions using Bioware software (Kistler) during a 60 s recording session at 500 Hz. Each animal was tested at least 3 times, with an interval of at least 30 min between trials. For each trial, we extracted the forces in the X-, Y- and Z-direction of the best push and calculated the highest peak resultant force. We used the highest peak resultant force across all trials for an animal as an estimate of its maximal push force.
External measurements were collected on all specimens used in the burrowing trials to characterize their external morphology. Head length from the tip of the snout to the back of the parietal bone, head width at the widest point of the head, head height at the tallest point of the head, lower jaw length measured from the tip of the jaw to the back of the retroarticular process, and body width at mid-body were measured using a digital calliper (Mitutoyo, ±0.1 mm). The snout–vent length was measured by stretching the animals along a ruler (±1 mm). Individuals were also weighed using an electronic balance or a spring scale (Ohaus, ±0.1 g) (Table 2).
Because species are not independent data points, their phylogeny was taken into account in our comparative analyses (Felsenstein, 1985). The phylogenetic tree of Jetz and Pyron (2018) was pruned to only include the species used in our study. Using 10,000 trees from VertLife.org, the maximum credibility tree was computed using the ‘maxCladeCred’ function in the Phangorn package in R (https://CRAN.R-project.org/package=phangorn).
3D geometric morphometrics
In addition to anatomical landmarks, 3D sliding semi-landmarks on curves were also used (Bardua et al., 2019a; Bookstein, 1991; Buser et al., 2018; Fabre et al., 2015, 2018; Gunz et al., 2005). Anatomical landmarks and sliding semi-landmarks were placed manually onto the crania and mandibles. All landmarks were placed by the same person (A.L.) using Stratovan Checkpoint (Stratovan corporation, v.2020.10.13.0859). Nineteen homologous landmarks and two curves were placed on the left mandible and 88 homologous landmarks and four curves were placed on the cranium (Figs 2 and 3 and Tables S2 and S3). Curves were resampled (see Botton-Divet et al., 2016, for a detailed description of the method) and all the sliding semi-landmarks were slid while minimizing the bending energy using the ‘slider3d’ function from the Morpho R package v.2.8 (https://CRAN.R-project.org/package=Morpho). Finally, a generalized Procrustes analysis (GPA) was performed using the ‘gpagen’ function from the Geomorph R package v.3.3.1 (https://CRAN.R-project.org/package=geomorph).
All the statistical analyses were performed in R version 4.0.3 (http://www.R-project.org/). The significance threshold was set at α=0.05. External measurements and forces were transformed logarithmically (log10) to fulfil assumptions of normality and homoscedasticity.
To assess the impact of size on shape, we performed a Procrustes regression on the GPA coordinates using the ‘procD.lm’ function from the Geomorph package. The log10 centroid size was used as a proxy for size. Residuals for both cranium and mandible were then computed and further referred to as allometry-free shapes, in order to examine shape variation not attributable to allometry.
To estimate the degree of similarity due to shared ancestry, a multivariate K-statistic (Adams, 2014) was calculated on the mean Procrustes coordinates of the cranium, the mandible and the external measurements using the ‘physignal’ function in the Geomorph package. A K-statistic was calculated on the forces using the ‘phylosig’ function in the Phytools package (https://CRAN.R-project.org/package=phytools). The phylogenetic signal was calculated under the assumption of Brownian motion (Blomberg et al., 2003). The higher the K-value, the stronger the phylogenetic signal. Values of K>1.0 describe data with a greater phylogenetic signal than expected, meaning that traits are conserved within the phylogeny.
To visualize the evolutionary patterns of shape variation in the cranium and mandible, we performed a principal component analysis (PCA) on the mean of the allometry-corrected shapes for each species using the ‘gm.prcomp’ function from the Geomorph package. Then, we projected the phylogeny onto the morphospace.
We then tested whether the force produced by secondarily aquatic typhlonectids and phylogenetically basal rhinatrematids differed from those of the other species, using an analysis of covariance (ANCOVA). As the force used by caecilians to push their heads into the soil is mainly generated by the body muscles (Bemis et al., 1983; Gaymer, 1971; Nussbaum and Naylor, 1982; O'Reilly et al., 1997), body width was used as a co-variate.
Next, a phylogenetic two-block partial least squares was performed using the ‘phylo.integration’ function from the Geomorph package with the resultant force as one block and the absolute external measurements as the second block to assess the co-variation between external head shape and force. As the analyses of the burrowing trials showed that forces are not exclusively generated in the antero-posterior direction (Table 2), we also explored the covariation between the absolute force vectors (Fx, Fy and Fz; block 1) and the cranium or mandible shape (block 2) using phylogenetic two-block partial least squares. Finally, phylogenetic generalized least squares (PGLS) regressions were performed in order to assess the relationship between maximum resultant push force and the shape of the cranium and mandible.
None of the measurements described in this paper (force measurements or external morphometrics) are considered procedures by French or European law. As such no ethics approval was required. All wild-caught animals were maintained for one night, checked for signs of stress and released at their exact site of capture. Captive animals were observed for signs of stress after experiments and monitored for signs of weight loss during the following week. None of the animals were harmed, or showed any signs of stress or weight loss after measurements.
The Procrustes regression of shape on log-transformed centroid size reveals that 12% of cranial shape variation is associated with size variation (P=0.001). An increase in cranium size across species is associated with the narrow and bullet-shaped cranium becoming wider and more triangular (Fig. S1). The impact of size on the mandibular shape is lower, as only 5% of the mandible shape variation is due to allometry (P=0.003). An increase in mandible size is associated with the mandible becoming thicker but with a shorter retroarticular process (Fig. S2).
The multivariate K-statistic calculated for the cranium (Kmult=0.83, P=0.001) and mandible (Kmult=0.88, P=0.001) was significant but the signal was only moderate and lower than expected by Brownian motion (K=1). However, no significant phylogenetic signal was detected for the external measurements (Kmult=0.7, P=0.39) or the resultant force (K=0.66, P=0.53).
Principal axes of shape variation of the cranium
The cranial morphological space described by the first five principal components (PCs) explains 73% of the total allometry-corrected shape variation. Each of the subsequent PCs explains less than 5% of the variation. PC1, explaining 23.3% of the total variation, mainly explains variation in the positioning of the mouth, with positive values corresponding to a terminal mouth and negative values to a subterminal mouth. A subterminal mouth is associated with a shortening of the maxillopalatine and the basal bones, resulting in a rounder and stockier cranium (compared with the more triangular and elongated shapes at the extreme positive PC1 scores). The braincase is also narrower towards the positive extreme of PC1, leading to a more open temporal region (zygokrotaphic skull), while the temporal region is closed on the negative extreme of PC1 (stegokrotaphic skull) (Fig. 4). PC2 (21.7% of the total variation) describes variation in the quadrate–squamosal complex. At the negative extreme lie species with narrow and posteriorly elongated quadrate and squamosal bones and an increased fenestration of the temporal region. Thus, species with zygokrotaphic skulls are associated mainly with negative PC2 scores, while species with stegokrotaphic skulls correspond to positive PC2 scores. Additionally, the overall cranium is flatter and more triangular towards the negative extreme of PC2, and rounder and more ellipsoidal towards the positive extreme (Fig. 4).
Principal axes of shape variation of the mandible
The mandibular morphological space on the first four PCs explains 80% of the total allometry-corrected shape variation. Each of the subsequent PCs explains less than 5% of the variation. PC1 (42.4% of the total variation) mainly corresponds to variation in the retroarticular process with positive scores corresponding to a retroarticular process lying in line with the pseudodentary. On the negative extreme of PC1, species show a dorso-medially curved retroarticular process, where the bending of the process is associated with a posterior shift of the articular surface. Additionally, the overall shape of the mandible is slender and elongate on the negative extreme. Bulkier mandible shapes are found at the positive extreme (Fig. 5). PC2 (19.2%) describes the overall curvature of the mandible. At the positive extreme, species show a straight pseudodentary and a retroarticular process almost in line with it. On the negative side, a strongly medially curved mandible is found (Fig. 5).
Variation in push force
Maximal push forces across the 13 species included in the study ranged from 0.53±0.16 N for Boulengerula fischeri, the smallest species in our dataset, to 16.11±4.96 N for Dermophis mexicanus, one of the biggest species in our dataset (see Table 2). The ANCOVA detected a significant effect of body width (F1,9=19.52; P=0.017) and group (F1,9=4.37; P=0.047) on maximal push force, suggesting that species from different groups (Typhlonectes, Rhinatrema, others) differed in their maximal push force, irrespective of variation in body width. Post hoc pairwise tests indicated that whereas the push force of Rhinatrema was not different from that in other caecilians (F1,9=17.78; P=0.10), individuals of the aquatic T. compressicauda produced significantly lower push forces than the other species (F1,9=5.99; P=0.037; Fig. 6).
Relationships between shape and force
The phylogenetic two-block partial least squares analyses (2B-PLS) show a significant and positive co-variation between external measurements and the resultant force (r-PLS=0.76; P=0.007). Head width, head height and body width drive most of the covariation with the resultant force (Fig. 7). There was no significant covariation between cranium shape (r-PLS=0.85; P=0.32) or mandible shape (r-PLS=0.71; P=0.75) and the different components of the push force (Fx, Fy and Fz).
Our PGLS regressions showed no significant relationship between shape and force, shape and size, or the interaction between force and size (see Table 3). Similarly, mandibular shape did not explain variation in resultant push force. However, mandibular shape was impacted by variation in body width (Table 3).
The objectives of this study were to: (1) assess the variability in skull shape across caecilians; (2) assess the variability in in vivo burrowing forces in caecilians; and (3) investigate the relationships between head shape, skull shape and in vivo burrowing forces.
Skull shape variability
Within the last decade, several studies have quantified the cranial shape variability among caecilians. These studies have shown that variation in cranial morphology is phylogenetically structured, yet cranial shape variability is rather high in caecilians (Bardua et al., 2019b; Sherratt et al., 2014). Our results confirm these observations. The main axes of variability are represented by features that are presumed to be related to burrowing performance: the position of the mouth, the temporal fenestration of the skull and the closure of orbits. Although these features do not always co-occur, there is a convergence towards a cranial morphology that has features that would be useful for efficient burrowing: closed orbits, stegokrotaphic skull, subterminal mouth and a globally narrow bullet-shaped skull. Interestingly, we also observed that the shape of the quadrate and the squamosal is highly variable across caecilians. These results are in line with observations by Bardua et al. (2019b) and Marshall et al. (2019) on the modularity of the caecilian skull. Among the 10-modules model proposed, the quadrate–squamosal module was the fastest evolving module (Bardua et al., 2019b). As suggested by Kleinteich et al. (2012), ventral bones are more likely to be shaped by demands on fossoriality rather than the bones on the dorsal surface of the skull. Here, we observed that the os basale is quite variable in length, and contributes to the curving of the ventral side of the skull. Although no 3D geometric morphometric studies have been done on caecilian mandibles to date, previous qualitative studies suggested that shape variation in the mandible is present mainly at the level of the retroarticular process (Wake, 2003). Our data confirm that the main axis of shape variation is the angle formed by the retroarticular process with the pseudodentary. This is interesting and suggests that the shape of the mandible is perhaps more constrained by feeding than burrowing given the important role of the curvature of the retroarticular process in generating bite force (Summers and Wake, 2005).
Some caecilians are known to use a combination of hydrostatic and internal concertina locomotion to move underground and to burrow (O'Reilly et al., 1997). These systems rely on the loose connection between the skin, associated with body wall muscles and the vertebral muscles (O'Reilly et al., 1997). As observed by Naylor and Nussbaum (1980) and Nussbaum (1977), rhinatrematids, the earliest diverging lineage of extant caecilians, possess a certain degree of independence between body wall and axial muscles, but less so compared with other species. This suggests that the ability to perform internal concertina was already present in the ancestor of caecilians. Our results show that the force exerted by Rhinatrema bivittatum is not significantly lower compared with other caecilians (Herrel and Measey, 2010; Summers and O'Reilly, 1997), suggesting that rhinatrematids are indeed capable of using internal concertina locomotion. Interestingly, however, their maximal push force for a given body diameter is somewhat lower than that of other terrestrial caecilians, suggesting that there may be subtle differences in the mechanics of burrowing.
In contrast, Typhlonectes is unable to perform internal concertina (Summers and O'Reilly, 1997). Typhlonectes natans were not able to move in tunnels equal to their body width and used normal concertina rather than internal concertina locomotion to move through wider tunnels (Summers and O'Reilly, 1997). Their X-ray images also suggest that axial and body wall muscles are not loosely connected, unlike in other terrestrial caecilians. Our results show that the push force produced by Typhlonectes compressicauda is significantly lower than that of other caecilians tested. Internal concertina locomotion thus appears to be a major component of force generation in caecilians (Ducey et al., 1993; O'Reilly et al., 1997). Typhlonectids thus probably lost the ability to perform internal concertina. Indeed, although aquatic caecilians are capable of head-first burrowing in soft substrates such as mud (Moodie, 1978; J.M. and M.W., personal observation), their aquatic locomotion relies on lateral undulation, for which a loose connection between body wall and axial muscles might be disadvantageous (Summers and O'Reilly, 1997).
However, our results do not support the predictions of Herrel and Measey (2010) suggesting that a higher degree of skin–vertebral independence should be correlated with higher burrowing forces. Indeed, our results show that species with a high degree of skin–vertebral independence, such as Herpele squalostoma and Schistometopum thomense, do not have greater burrowing forces compared with species with a lower degree of independence, such as Geotrypetes seraphini and Boulengerula taitanus. However, the skin vertebral independence observed by Summers and O'Reilly (1997) in Dermophis mexicanus, 5 times greater than the values observed by Herrel and Measey (2010) for Schistometopum thomense, corresponds to the greatest burrowing force measured in our dataset. Additional X-ray data in combination with direct measurements of forces are needed to fully understand the relationship between body wall and axial muscle independence and push force in caecilians.
Relationships between shape and burrowing force
Digging is an energetically demanding process (Wu et al., 2015). The costs of burrowing increase exponentially with increasing body diameter (Gans, 1968; Navas et al., 2004; Vanhooydonck et al., 2011). The most energetically efficient way to increase muscle force, and thus maximal push force, should be to increase body length rather than body diameter as suggested by previous studies (Barros et al., 2011; Gans, 1968; Herrel et al., 2021; Hohl et al., 2017; Le Guilloux et al., 2020; Measey and Herrel, 2006; Navas et al., 2004; Vanhooydonck et al., 2011). In our study, all the external measurements used to characterize external head shape are positively correlated with maximal in vivo burrowing force. Overall, bigger and specifically wider animals are capable of pushing harder (e.g. scolecophidians: Herrel et al., 2021; skinks: Le Guilloux et al., 2020; amphisbaenians: Navas et al., 2004). Indeed, among the variables tested, body width and head width are most strongly correlated with push force. These results are in line with previous studies on limbless burrowing tetrapods: in skinks and amphisbaenians, absolute push force is best predicted by head width and body diameter (Le Guilloux et al., 2020; Navas et al., 2004). As the muscles responsible for internal concertina and hydrostatic locomotion in caecilians are positioned to the lateral side of the body (Gaymer, 1971; O'Reilly et al., 1997), it makes sense that body diameter is positively correlated with maximal push force and probably reflects the presence of muscles with larger cross-sectional areas. This further suggests that burrowing energetics are probably not an important driver of overall shape.
Although the cranial shape variation in caecilians is mainly explained by features thought to be related to an increased burrowing performance, our results show that there is no relationship between skull shape and maximal in vivo burrowing force. These results are in line with the predictions of Kleinteich et al. (2012) based on biomechanical models: complete coverage of the temporal region (stegokrotaphy) is not associated with an increase in performance during burrowing. As observed both by Kleinteich et al. (2012) and in our study, the os basale is one of the drivers of shape variation in caecilian skulls. This variation in length of the os basale is associated with an overall curving of the ventral surface of the skull. As suggested by Kleinteich et al. (2012), there is an optimal head angle for burrowing around which the strain and stress are minimal. Although we found no evident relationship between maximal burrowing force and shape variation, other parameters of burrowing such as head angle or burrowing speed might reveal a positive correlation with skull shape and need to be explored further.
In addition to a role in head-first burrowing, the cranium has multiple functions such as feeding and gas exchange, and houses the brain and major sensory organs. The skull is thus expected to be shaped by many competing functional demands. As suggested by previous authors, skull shape might be related to selective pressures associated with feeding rather than burrowing. The presence of a mobile quadrate (streptostyly) has been suggested to impact bite force generation (Kleinteich et al., 2008b; Summers and Wake, 2005; Wake and Hanken, 1982). Moreover, variation in caecilian skulls appears to be concentrated in the temporal region (Bardua et al., 2019b). As the fenestration of the temporal region does not seem to be related to burrowing performance (Kleinteich et al., 2012), this region may be constrained by the insertion of the jaw adductor muscles. Caecilians possess a unique dual jaw-closing mechanism (Nussbaum, 1983) involving the m. adductor mandibulae and the m. interhyoideus posterior. Although the latter inserts onto the retroarticular process, the m. adductor mandibulae takes its origin from the dorsolateral surface of the parietal bone and inserts on the mandible. A complete coverage of the temporal region might then play a role in providing additional muscle insertion area. However, a complete coverage of the temporal region could also restrain the space available for the adductors.
Moreover, although the mandibles are thought to play a preponderant role in burrowing (terminal versus sub-terminal mouth), we found no relationship between mandible shape and maximal in vivo burrowing force in caecilians. However, as most of the shape variation occurs around the retroarticular process, mandible shape will probably have consequences for the biomechanics of jaw closure (Summers and Wake, 2005; Wake, 2003). More data on burrowing and feeding performance, such as burrowing speed and bite force or long-axis twisting during the reduction of large prey (Measey and Herrel, 2006), are needed to better understand the selective pressures that have shaped the skull in caecilians.
Our results suggest that differences in the degree of fossoriality are not associated with the evolution of skull shape in caecilians. Moreover, our results confirm previous studies suggesting that: (1) independence between body wall and axial muscles is necessary to perform internal concertina locomotion, and (2) this locomotor mode is the main agent of force production in caecilians. However, more performance data are needed to fully understand the selective pressures that have shaped the skull in caecilians.
We thank I. Josipovic and the people at Centre for X-Ray Tomography at Ghent University for their help with CT scanning. We thank the Museum of Zoology (University of Michigan), the Amphibian & Reptile Diversity Research Center (University of Texas Arlington), the Royal Museum of Central Africa (Brussels) and the Zoological Museum (Hamburg) and their curators for the loan of some key specimens. We also thank MorphoSource for making scans available. We would also thank Alexander Kupfer and an anonymous reviewer for helpful and constructive comments on a previous version of the manuscript. A.L. thanks A.-C. Fabre and M. Zelditch for their helpful discussions on the statistical analyses. A.L. also thanks M. H. Wake for the gift of Dermophis specimens and her helpful discussion on the manuscript.
Conceptualization: A.L., A.H., D.A.; Methodology: A.L.; Formal analysis: A.L.; Investigation: A.L., A.H.; Resources: A.L., B.D., M.W., J.M., J.C.O., N.J.K., P.G., J.B., T.K., A.H.; Data curation: A.L., B.D.; Writing - original draft: A.L.; Writing - review & editing: A.L., B.D., M.W., J.M., J.C.O., N.J.K., P.G., J.B., T.K., L.V., A.H., D.A.; Visualization: A.L., D.A.; Supervision: A.H., D.A.; Project administration: A.H., D.A.; Funding acquisition: A.L., J.M., L.V.
This study was supported by the Research Foundation, Flanders (Fonds Wetenschappelijk Onderzoek, grant 11D5819N), a Tournesol travel grant and a European Union Marie Curie Fellowship (HPMF-CT-2001-01407), field work and visiting fellowship of the Fonds Wetenschappelijk Onderzoek, Flanders, Belgium (FWO-Vl) to J.M. The special research fund of Ghent University (BOF-UGent) is acknowledged for financial support of the UGCT Centre of Expertise (BOF.EXP.2017.0007).
The authors declare no competing or financial interests.