Is solid always best? Cranial performance in solid and fenestrated caecilian skulls

Caecilians are a monophyletic group of amphibians with skulls that are thought to be highly specialized for their burrowing lifestyle; they are wedge shaped, compact and robust (Wake, 1993). Many bones of the skull have been fused into larger compound elements, e.g. the os basale of caecilians comprises all of the bones of the skull base [i.e. the parasphenoid, the exoccipitals, and the caudal parts of the neurocranium, including the otic capsules (Wake, 2003) (Fig. 1)]. There are two distinct skull types in caecilians: (1) zygokrotaphic, in which the skull is fenestrated between the squamosal and the parietal, and (2) stegokrotaphic, in which the skull is completely roofed.

Based on recently published caecilian phylogenies (Roelants et al., 2007; Zhang and Wake, 2009; Pyron and Wiens, 2011; Wilkinson et al., 2011) (Fig. 2), the zygokrotaphic skull has evolved independently several times in caecilians, in the Scolecomorphidae, Typhlonectidae and Dermophiidae (Brand, 1956; Taylor, 1969; Nussbaum, 1977; Nussbaum, 1985; Wilkinson and Nussbaum, 1997; Müller et al., 2009). Zygokrotaphy in the Rhinatrematidae, however, has usually been considered to be the ancestral condition for the Gymnophiona (Nussbaum, 1977; Nussbaum, 1983; Wake, 2003; Müller, 2007), and the reduction of bone coverage in the temporal skull region was considered to be homologous among caecilians, frogs and salamanders. Thus, the completely roofed stegokrotaphic skulls of some caecilians were considered to be secondarily derived in association with their burrowing lifestyle (Peter, 1898; Goodrich, 1958; Parsons and Williams, 1963; Nussbaum, 1983; but see Marcus et al., 1933; Carroll and Currie, 1975). This view is supported by the absence of a real suture between the squamosal and the parietal in the stegokrotaphic skull. However, the discovery of a stem group caecilian with a stegokrotaphic skull, Eocaecilia micropodia, has challenged this view and suggests that the zygokrotaphic skull in the Rhinatrematidae is a derived condition and that stegokrotaphy is ancestral for all caecilians (Jenkins and Walsh, 1993; Carroll, 2000; Jenkins et al., 2007).

Despite coverage of the temporal region by the squamosal in the stegokrotaphic caecilian skull type, there is always a narrow gap between the squamosal and the parietal instead of a tight suture (Wiedersheim, 1879; Sarasin and Sarasin, 1887–1890; Peter, 1898; Versluys, 1912; Abel, 1919; Marcus et al., 1933; Goodrich, 1958; Lawson, 1963; Taylor, 1969; Nussbaum, 1977; Nussbaum, 1983; Wake and Hanken, 1982). This gap is supposed to allow movement of the squamosal and the attached quadrate and thus plays a role in a complex cranial kinesis (Versluys, 1912; Marcus et al., 1933; Iordansky, 1990; Iordansky, 2000). Movement of the squamosal is related to the unique caecilian dual jaw closing mechanism that involves an accessory ventral jaw closing muscle (m. interhyoideus posterior) acting simultaneously with the primary jaw closing muscles (mm. levatores mandibulae) (Bemis et al., 1983; Nussbaum, 1983). The movement of the squamosal during feeding, as revealed by three-dimensional (3-D) modeling, is a small mediolateral rocking that does not expose a substantial gap in the skull (Kleinteich et al., 2008; Kleinteich, 2010). The stem group caecilian E. micropodia, however, differs from extant caecilians by having additional bones in the temporal region (e.g. a presumed tabular) and a different jaw joint that is flat instead of a deep groove (Jenkins et al., 2007), which suggests the presence of a different jaw closing mechanism, possibly without movements of the squamosal.

Understanding the costs and benefits of one skull type over the other is crucial for evaluating amphibian skull evolution because it gives insight into how strongly the biology of a caecilian is affected by its skull. The stegokrotaphic skull would, on casual inspection, appear better suited to burrowing than the zygokrotaphic skull with the large unroofed region and, for two caecilian groups with zygokrotaphic skulls, there might actually be an ecological link to poorer burrowing performance; rhinatrematids are considered to be poor burrowers that can be readily trapped during surface activity (Nussbaum, 1983; Gower et al., 2010); typhlonectids are mainly aquatic species and may therefore be relieved of the constraints of burrowing (Wilkinson and Nussbaum, 1997). However, the zygokrotaphic scolecomorphids are considered to be specialized burrowers (Nussbaum, 1977; Nussbaum, 1983); furthermore, the zygokrotaphic dermophiid caecilian Geotrypetes seraphini shows similar burrowing performance to other stegokrotaphic caecilians (Herrel and Measey, 2010).

Here, we present an experimental study where we modified the 3-D geometry of zygokrotaphic and stegokrotaphic caecilian skulls to artificially render zygokrotaphic skulls stegokrotaphic and vice versa. We further modified the stegokrotaphic skull models to test a third, akinetic skull type. The original and modified skull geometries were tested using finite element analysis (FEA) for their performance under the frontal loading regime that caecilians encounter during burrowing. The skulls were oriented at different angles relative to their anterior–posterior axis to simulate varying directions of the frontal loads that caecilians are likely to encounter by dorso-ventral movements of the head to manipulate the substrate (Wake, 1993). The aims of this study were: (1) to evaluate the sensitivity of the strain distribution in caecilian skulls as the head moves through different angles during burrowing; (2) to quantify the difference in deformation of zygokrotaphic and stegokrotaphic caecilian skulls under frontal loading; (3) to determine whether deflection during burrowing is affected by the state of the skull roofing by generating and testing hypothetical morphologies in which the skull roofing state is modified relative to the condition present in the actual species; and (4) to test whether deformation during burrowing in an akinetic skull representative of a putative caecilian ancestor is affected by the restricted movements of the squamosal.

Six adult caecilian specimens from six different species were available for this study (Table 1). The specimens were made available by the Zoological Museum Hamburg, Germany (ZMH), Mark Wilkinson (MW; Natural History Museum London, UK), Alexander Kupfer (AK; University of Siegen, Germany), and the Muséum National d’Histoire Naturelle, Paris, France (MNHN). Our sampling comprised three species with zygokrotaphic skulls and three species with stegokrotaphic skulls (Table 1). Further, the species were sampled broadly across caecilian phylogeny and represent six of the nine families that are recognized by Wilkinson et al. (Wilkinson et al., 2011) (Fig. 2).

The species Ichthyophis kohtaoensis was originally described as endemic to Koh Tao Island in Thailand (Taylor, 1960). It is currently unclear if I. kohtaoensis can also be found on the mainland of Thailand or if the mainland populations belong to a different species. The specimen examined herein was derived from the pet trade and it is not possible to assign it to either the island or the mainland populations of the genus Ichthyophis of Thailand with certainty and thus we prefer to refer to this specimen as Ichthyophis cf. kohtaoensis.

We obtained 3-D volume datasets with high-resolution micro-computed tomography (HRμCT). Except for the Rhinatrema bivittatum specimen (see below), HRμCT scanning was performed with synchrotron-based x-ray radiation at the beamline W2 of the DORIS III accelerator ring of the German Electron Synchrotron (DESY) in Hamburg, Germany. This beamline is operated by the Helmholtz Center Geesthacht, Germany. For synchrotron-based x-ray radiation HRμCT, we used monochromatic x-rays and adjusted the energy of the x-ray beam individually to the samples (Table 1). Details of the setup for HRμCT imaging at beamline W2 of the German Electron Synchrotron have been published previously (Beckmann et al., 2006; Kleinteich et al., 2008). The R. bivittatum specimen was HRμCT-scanned at the 3-D Morphometrics Laboratory facility at the University of Calgary, Canada. The scan was performed on a Scanco μCT35 scanner (Scanco Medical AG, Brüttisellen, Switzerland) with the beam energy set at 55 kVp and 72 μA, and a voxel size of 20 μm.

The Boulengerula taitana and R. bivittatum specimens were HRμCT scanned in 70% ethanol; the remainder specimens were decapitated, the heads were exposed to –80°C for 3 h and then freeze dried with a Lyovac GT2 freeze drying system (Leybold-Heraeus GmbH, Hanau, Germany) for 24 h according to the procedure described previously (Meryman, 1960; Meryman, 1961) prior to HRμCT imaging. Freeze drying is known to increase the contrast in x-ray images (Follett, 1968). Compared with other drying methods, freeze-drying results in less volume shrinkage (Boyde, 1978) and, especially for hard tissues, shrinkage effects due to freeze drying are considered to be negligible.

HRμCT imaging at DESY resulted in volume datasets with resolutions from 3.9 to 9.2 μm; the dataset of the R. bivittatum specimen had a resolution of 20.0 μm (Table 1). Because of limitations in the computer hardware for visualization of synchrotron HRμCT datasets, we pairwise combined neighboring voxels in all three dimensions of the original synchrotron HTμCT datasets (two-fold binning), which decreased the resolution of the datasets by a factor of two.

The HRμCT datasets were then loaded into Amira 5.2 (Visage Imaging GmbH, Berlin, Germany). We segmented the cranium except for the lower jaw and the stapes from the volumetric dataset by using the LabelField function in Amira. Segmentation of HRμCT data with the LabelField function results in a so-called Labels dataset that contains the assignments of voxels in the HRμCT data to specified Materials, in our case the only Material in the Labels dataset was the segmented cranium. Based on the Labels dataset derived from the original anatomy, we generated two artificial anatomies for each specimen: (1) by removing voxels in the Labels dataset from the dorsal region of the squamosal in stegokrotaphic species (i.e. converting them to a zygokrotaphic skull type) or by adding voxels to the dorsal edge of the squamosal in zygokrotaphic species (i.e. converting them to stegokrotaphic skulls), and (2) by filling the gap between the squamosal and the parietal with voxels to make a tight connection of the squamosal to the remainder cranium (i.e. converting them to akinetic skulls) (Fig. 3). Other than in previously published studies on manipulations of 3-D datasets (Stayton, 2009; Stayton, 2011), our approach is not based on landmark data and the geometric morphometric method.

The Labels datasets were then used to calculate surface models of the skulls. Each species was represented by three different surface models based on the original and modified Labels datasets (Fig. 3). Surfaces that are generated with Amira contain an unnecessarily high number of polygons, and a reduction in the number of polygons usually has no or only little effect on the appearance of the surface model. We reduced the number of polygons until we could visually detect a decrease in the quality of the surface representation of our HRμCT datasets. Surfaces were edited in Amira to remove triangular geometries with low aspect ratios and intersecting polygons. We scaled all surface models to identical surface areas to allow for direct comparisons of calculated stresses between the different skulls (Dumont et al., 2009). Original surface areas and linear scaling factors that we applied in the x, y and z directions are shown in supplementary material Table S1.

FEA was performed with Marc Mentat 2005 R3 (MSC Software Corporation, Santa Ana, CA, USA) and the surface models were imported with a 3D solid geometry from Amira. The models consisted of approximately 2.3 million to 3.2 million tetrahedral elements (supplementary material Table S1). We treated the elements in the models as an isotropic material with a Young’s modulus of 10 GPa and a Poisson’s ratio of 0.3. These values have been used in previous studies that involved FEA to assess patterns in non-mammalian tetrapod skull evolution (e.g. Moazen et al., 2008; Moazen et al., 2009). Further, Dumont et al. (Dumont et al., 2009) demonstrated that FEA can be a useful tool for comparisons of shapes, even in the absence of knowledge on the exact material attributes in vivo.

We applied two constraints (in Marc Mentat Boundary Conditions) to the models (Fig. 1): (1) the nodes that define the joint area of the occipital condyles were prevented from movements by applying a Fixed Displacement of zero for translations and rotations along each axis to them (Fig. 1, green areas); (2) nodes in the area around the nasal capsule that define the rostral surface of the skull were selected and we applied Point Loads to them that summed up to a force of 10 N (Fig. 1, pink areas). A burrowing force of 10 N lies well within the range of forward pushing forces that were reported for caecilians; in vivo measurements showed that peak forward pushing forces of caecilians can be up to 20 N (O’Reilly et al., 1997). supplementary material Table S1 shows the numbers of constrained nodes; the total number of constrained nodes is less than 1% of the available nodes in the models, and small variations in the number of constrained nodes between different models are assumed to be negligible. To account for different angles of the head relative to the axis along which caecilians are moving forward, we incrementally rotated the skulls in five-degree steps for –15 deg, –10 deg, –5 deg, 0 deg, 5 deg, 10 deg, 15 deg relative to the neutral axis of the skull (as defined by the occipital condyles and the rostral tip of the nasal capsule).

We used three different measurements to describe the performance of the different skull models under the loading regime we applied herein: (1) total strain energy, (2) maximum Von Mises stress and (3) mean Von Mises stress. Total strain energy equals the sum of all deformations of the single elements in a model, multiplied by the applied force. This measurement reflects the amount of work that is put into elastic deformation of the skull. We can assume that the energy that caecilians invest to load their skulls is meant to push the skulls forward into the substrate and not to store elastic energy by bending the skulls. Thus, caecilians are expected to minimize total strain energy on their skulls, i.e. to be more energy efficient during burrowing.

Maximum Von Mises stress is a performance measurement that indicates how close the models are to failure. Structures that encounter higher Von Mises stress are closer to failure than structures with lower values for Von Mises stress. However, it is known that the constraints on finite element models (i.e. Fixed Displacements and Point Loads) can cause artificially high maximum Von Mises stresses in the regions of the constraints, which can confound the interpretation of maximum Von Mises stress values (Dumont et al., 2009). Further, because maximum Von Mises stress only represents one node in the model, the key value can refer to different regions in different specimens or even different regions in the same specimen under different loading conditions.

We calculated mean Von Mises stresses under the assumption that skulls that perform better during burrowing will show a decrease in Von Mises stress over all elements. Mean Von Mises stress is supposed to be less affected by the constraints of the finite element models and accounts for all nodes in the model.

While Von Mises stress scales linearly with the force per surface area ratio in finite element models, total strain energy scales with the force per volume ratio (Dumont et al., 2009). Although the finite element models herein had similar force per surface area ratios, their volumes differed slightly due to different cranial shapes. We calculated a force scale factor (supplementary material Table S1) that is based on the differences in the volume of the models. We then multiplied the total strain energies that resulted from the FEA with the squared force scale factors [because total strain energy scales with the square of the force (Dumont et al., 2009)] to calculate total strain energies for equally scaled force per volume ratios.

Total strain energy and Von Mises stresses are a standard output for FEA with Marc Mentat. Values for Von Mises stresses for every node were imported to the statistical computing environment R 2.13.1 (http://www.r-project.org) to calculate mean Von Mises stress. We fitted the distribution of total strain energy and mean Von Mises stress to a raw polynomial function of second order over the different head angles by using a linear model in R. From the resulting regression equation, we calculated the angles at which both performance measurements were minimal, as well as the angles where total strain energy and mean Von Mises stress were increased by 10, 20 and 30% compared with the minimum. For maximum Von Mises stress we did not calculate a polynomial regression because values for maximum Von Mises stress at different head angles actually refer to different nodes and thus to different regions in the models.

The effects of skull type and head angle on total strain energy, maximum Von Mises stress and mean Von Mises stress were evaluated by calculating a two-way analysis of variance (ANOVA) for each performance measurement. Further, we grouped stegokrotaphic and zygokrotaphic skulls as kinetic skulls and compared them with akinetic skulls by using three-way ANOVAs that considered total strain energy, mean Von Mises stress and maximum Von Mises stress as dependent on the presence of cranial kinesis, the species and the head angle. All statistical calculations were performed within the statistical computing environment R 2.13.1.

Application of a frontal load that is parallel to the neutral axis of the skull of caecilians (i.e. at a 0 deg head angle) causes a slight posteroventral bending of the nasal capsule region (Fig. 4). Variation in the head angle will alter this pattern and either increase this effect if the nose is lowered relative to the occipital condyles (i.e. positive head angles herein) or cause the nasal capsule region to be pushed dorsally instead of ventrally if the nose is elevated compared with the occipital condyles (i.e. negative head angles herein).

Starting from the neutral axis (0 deg), mean Von Mises stress (Fig. 5, supplementary material Table S2) and total strain energy (Fig. 6, supplementary material Table S3) will increase as either the nose is lifted relative to the occipital condyles or vice versa. Mean Von Mises stress and total strain energy differ significantly at different head angles (mean Von Mises stress, d.f.=6, F=318.23, P<0.0001; total strain energy, d.f.=6, F=152.96, P<0.0001). A linear model fit of the observed distributions of mean Von Mises stress (Fig. 5) and total strain energy (Fig. 6) over the different head angles onto a polynomial function of the second order results in high correlation coefficients (R²>0.97) and a significant relationship between mean Von Mises stress and head angle and between total strain energy and head angle (P<0.00045) for all models tested herein (supplementary material Table S4). This shows that the correlation between head angle and skull performance during burrowing can be described as a parabola.

Fig. 7 shows the distributions of Von Mises stress over the skulls of caecilians in dorsal and ventral view at two different head angles relative to the neutral axis of the skulls. At –10 deg (i.e. the nose is elevated relative to the occipital condyles) and +10 deg (i.e. the nose is lowered relative to the occipital condyles), the highest Von Mises stresses are concentrated in the caudal parts of the skulls. High stresses occur around the occipital condyles and the caudal parts of the os basale that comprise the otic capsules. At +10 deg, Von Mises stresses are much higher on the ventral surface of the skulls than at –10 deg. The area of high Von Mises stress that is concentrated around the premaxillaries and the vomers at –10 deg extends all over the ventral surface of the skull from the premaxillary to the caudal parts of the os basale at +10 deg. Further, compared with the situation at –10 deg, the dorsal skull bones encounter higher Von Mises stresses at +10 deg, except for in Typhlonectes natans (Fig. 7).

Based on the results of the polynomial regression, we calculated the head angles that were optimal for burrowing as well as the angles for which mean Von Mises stress and total strain energy increased by 10, 20 and 30% (Table 2). In all models that we tested, the optimal head angle for burrowing was slightly off the neutral axis in a negative direction, i.e. with a slightly elevated nose. Optimal head angles to minimize mean Von Mises stress were calculated to lie between –6.1 deg (stegokrotaphic R. bivittatum) and –2.5 deg (stegokrotaphic and zygokrotaphic T. natans); for total strain energy, the optimal angles were found to range from –4.1 deg (stegokrotaphic I. kohtaoensis) to –1.7 deg (akinetic B. taitana) (Table 2). The offset angles from the optimum that result in an increase of mean Von Mises stress by 30% were calculated to lie between 7.8 deg and 9.2 deg; for a 30% increase in total strain energy the offset angles were even smaller and ranged from 2.5 deg to 3.8 deg (Table 2). In both performance parameters, B. taitana and T. natans show the smallest offset angles.

Fig. 8 shows the distribution of Von Mises stress on all the different models we tested herein at the optimal angle (based on mean Von Mises stress; Table 2, Fig. 5). Comparable to the situation at –10 deg and +10 deg head angle, the highest Von Mises stresses occur at the premaxillaries and vomers, and close to the occipital condyles. However, at the optimal angle, Von Mises stress in the caudal parts and on the dorsal surface of the skull is lower compared with at –10 deg and +10 deg (Fig. 7). Von Mises stress on the ventral surface of the skulls at the optimal angle (Fig. 8) is higher compared with those at a head angle of –10 deg (Fig. 7) and extends over a wider area; compared with the situation at+10 deg (Fig. 7), the ventral surface of the skulls experiences lower Von Mises stresses at the optimal angle (Fig. 8).

For maximum Von Mises stress, the observed values actually refer to different nodes in different regions of the same model. Thus, the collected data on maximum Von Mises stress at different head angles of the same model is not strictly continuous and we did not calculate regression equations for maximum Von Mises stress. However, Fig. 9 and supplementary material Table S5 show the values for maximum Von Mises stress for the models we tested herein at varying head angles. Comparable to mean Von Mises stress and total strain energy, maximum Von Mises stress tends to increase if the nose is lifted beyond an optimal angle that appears to lie between –10 deg and 0 deg and also increases if the nose is lowered above the optimal angle. The three tested models that were derived from the R. bivittatum HRμCT dataset, however, show a steady decrease in maximum Von Mises stress with increasing head angle, which is likely to be an artifact that is caused by a high point load around a few nodes in these models. Exclusion of these nodes results in the same pattern that we observed in the other caecilian species.

The three skull types showed very similar behaviors under the frontal loading regime that we applied herein. Values for total strain energy (d.f.=2, F=0.51, P=0.60), maximum Von Mises stress (d.f.=2, F=0.37, P=0.69) and mean Von Mises stress (d.f.=2, F=0.49, P=0.62) were all found to not be significantly different for the three skull types tested.

Fig. 8 compares the distributions of Von Mises stress on the dorsal and the ventral surfaces of the different skull models for each species with each other. In none of the six caecilian species do the dorsal parts of the squamosal experience notable Von Mises stresses under a frontal loading regime. Thus, varying the shape of the squamosal by removing or adding material does not have an effect on the distribution of Von Mises stress over the skull. Further, because the squamosal is not under stress, the shape of the squamosal, and thus the skull type, has no effect on the response of the skull to different head angles (Figs 5, 6, 9). The interactions of skull type with head angle were found to be not statistically significant for total strain energy (d.f.=12, F=0.03, P>0.99), maximum Von Mises stress (d.f.=12, F=0.01, P>0.99), and mean Von Mises stress (d.f.=12, F=0.03, P>0.99).

However, akinetic skulls experience slightly higher Von Mises stresses at the tooth-bearing bones that define the outline of the skulls in ventral view (i.e. the lateral regions of the premaxillary and the maxillopalatine) and slightly lower Von Mises stresses at the central ventral bones compared with stegokrotaphic and zygokrotaphic skulls (Fig. 8). By combining zygokrotaphic and stegokrotaphic skulls as kinetic skull morphologies and comparing them to the akinetic skull models in a three-way ANOVA that accounts for presence of skull kinesis, head angle and species, we found that total strain energy was significantly different between kinetic and akinetic skulls (d.f.=1, F=27.6, P<0.0001). In all the species we investigated herein, we found consistently lower (by approximately 5–15%) values for total strain energy in the akinetic skulls compared with kinetic skulls (supplementary material Table S3). Mean Von Mises stress (d.f.=1, F=1.53, P=0.22) and maximum Von Mises stress (d.f.=1, F=0.01, P=0.91) were found to not significantly differ between kinetic and akinetic skulls.

Our results demonstrate that the two different caecilian skull architectures (i.e. stegokrotaphy and zygokrotaphy) do not have an effect on skull performance under a frontal loading regime, which caecilians encounter while pushing the skull into the substrate during burrowing. Further, the kinetic caecilian skulls show similar patterns of Von Mises stress distribution to akinetic skulls, which may have been present in Paleozoic caecilian and tetrapod ancestors. However, we found that akinetic skulls have consistently lower values for total strain energy. We further present evidence that there is an optimal head angle for burrowing in caecilians around which total strain energy and mean Von Mises stress are minimal. At the optimal head angle, stress close to the occipital condyles is minimized, presumably due to a reduction in torque around the occipital joint. We also show that most of the load that caecilian skulls encounter during burrowing is transmitted along the ventral surface of the skull. Thus, the ventral skull bones are more likely to be shaped by demands of fossoriality than the bones on the dorsal surface of the skull.

Traditionally, there has been an emphasis on the different caecilian skull types as they relate to the question of whether the solid skull roof (stegokrotaphy) is the ancestral condition (Marcus et al., 1933; Carroll and Currie, 1975) or whether it is derived as an adaptation to the burrowing lifestyle of caecilians (Peter, 1898; Abel, 1919; Goodrich, 1958; Parsons and Williams, 1963). There is no strong consensus among authors as to which of the two caecilian skull types is ancestral (Wake, 2003; Carroll, 2007; Carroll, 2009; Jenkins et al., 2007). Although we did not gather data that directly bears on the question of caecilian ancestry [e.g. temnospondyl versus lepospondyl origins (see Anderson, 2008)], we have shown that with respect to the current utility of the skull during burrowing, a stegokrotaphic skull does not result in improved performance over the zygokrotaphic skull. This suggests that if stegokrotaphy is the derived condition within caecilians, it did not evolve primarily as an adaptation to burrowing. Additionally, if stegokrotaphy is the ancestral condition, the evolution of temporal fenestration (zygokrotaphy) does not imply a coupled decrease in cranial performance during burrowing. No matter whether stegokrotaphy or zygokrotaphy is ancestral for caecilians, zygokrotaphy evolved at least three times independently, i.e. within the Scolecomorphidae, Typhlonectidae and Dermophiidae. If zygokrotaphy in the Rhinatrematidae is the derived condition, it may have even evolved four times (Fig. 2). The lack of an effect on the most constraining function of the skull, its use for burrowing, may have facilitated the multiple evolution of zygokrotaphy.

Zygokrotaphy very likely results in more space for the jaw closing muscles (mm. levatores mandibulae) and thus might be related to differences in feeding biomechanics rather than to differences in burrowing performance. In the zygokrotaphic skulls of the Rhinatrematidae, the jaw closing muscles pass through the temporal fenestrae and reach towards the dorsal midline of the skull roof; in all other caecilians with zygokrotaphic skulls, the jaw closing muscles do not extend beyond the dorsal edge of the squamosal (Nussbaum, 1977; Nussbaum, 1983; Wake, 2003). However, we would like to point out that although the jaw closing muscles in caecilians with secondarily derived zygokrotaphic skulls do not extend as far dorsal as in species of the Rhinatrematidae, zygokrotaphy does leave more space for the jaw closing muscles. Kleinteich et al. showed that the jaw closing muscles that act on the zygokrotaphic skull of T. natans have a higher contribution to total bite force than the jaw closing muscles in the stegokrotaphic species B. taitana, Siphonops annulatus and I. cf. Kohtaoensis (Kleinteich et al., 2008; Kleinteich, 2010). In the zygokrotaphic G. seraphini, the jaw closing muscles extend onto the dorsal surface of the parietal (T.K., personal observation), while in stegokrotaphic caecilians, the jaw closing muscles originate from the lateral edge of the parietal (Edgeworth, 1935; Lawson, 1965; Iordansky, 1996; Iordansky, 2010).

In all species that we tested, independent of the skull type, we found that Von Mises stress was higher on the ventral surface of the skull than on the dorsal surface. This is because at an optimal head angle for burrowing, the sequence of the ventral bones (i.e. premaxillary, vomer and parasphenoid portion of the os basale) is almost in line with an axis that directly connects the nasal capsule with the occipital condyles. High Von Mises stresses on the ventral surface of the skull are caused by a slight downward bending of the skull under load (Fig. 4), causing the ventral surface of the skull to be under compression. The ventral bending of the skull follows the pattern of cranial kinesis that was previously predicted for caecilians by Iordansky (Iordansky, 1990; Iordansky, 2000). However, although much attention has been paid to the skull roof as an adaptation for burrowing (Peter, 1898; Nussbaum, 1977; Wake and Hanken, 1982; Bemis et al., 1983; Nussbaum, 1983; Wake, 1993), the role of the palate during burrowing has been neglected. Preliminary examination of the nature of the sutures between palatal bones reveals broadly overlapping joints filled by a dense network of collagen fibers, suggesting at least minimal kinesis in this region (H.C.M., personal observation). The anatomy of these sutures is consistent with a compression-sink model that would alleviate stress caused under a compression scenario (Wake, 1993). This hypothesis has yet to be tested empirically and is currently under investigation.

Our results show that a skull that is akinetic but otherwise similar to a skull of extant caecilians has no major advantages or disadvantages when exposed to the conditions mimicking burrowing. However, the plots of Von Mises stress reveal slightly decreased Von Mises stress on the vomer and the os basale and a slightly increased Von Mises stress on the tooth-bearing bones. Further, akinetic skulls bend less, which is indicated by significantly lower total strain energy. There might be a trade-off between reduction of skull bending and loading of the bones that define the ventral outline of a caecilian skull. The laterally exposed parts of the premaxillary, the maxillopalatine and also the squamosal tend to be much thinner bones than the central bones, especially the os basale. We hypothesize that, besides the advantages of a kinetic skull during biting that were discussed in detail by Summers and Wake (Summers and Wake, 2005) and Kleinteich et al. (Kleinteich et al., 2008), detachment of the squamosal from the dorsal dermal bones also might be a mechanism to guide loads towards the ventral centre of the skulls where thicker bones, and extensive overlapping (and possibly weakly kinetic) joints, are present.

Stress on the ventral surface in kinetic and akinetic skulls can be reduced by raising the nose slightly compared with the occipital condyles. This reduces ventral bending and at some point will push the tip of the nose dorsally rather than ventrally. However, this results in an increase in Von Mises stress on the dorsal surface of the skull and in the caudal region of the os basale. Raising the nose beyond an optimal angle results in higher mean Von Mises stresses, higher maximum Von Mises stresses and higher total strain energies. On the other hand, lowering the nose beyond an optimal angle during burrowing results in higher Von Mises stresses on the ventral and dorsal surface of the skulls and at the occipital condyles and also decreases overall performance. Even slight deviations from the optimal head angle in either direction might be of biological significance. This leads us to predict that dorso-ventral movements of the skull during burrowing will be restricted to a narrow range near the zero angle of incidence. Although the concept of optimality in head angles during burrowing to reduce the stress and bending energies on the skulls is intuitive, it has not been tested in vivo. Future studies on skull loads during burrowing in vertebrates should address the sensitivity of the results to slight variations in the head angle.

We are grateful to Mark Wilkinson (Natural History Museum, London, UK) and Alexander Kupfer (Universität Siegen, Germany), who kindly donated specimens for this study. We appreciate the help of Tamer Fawzy (Hamburg, Germany), Thomas Dejaco (Innsbruck, Austria) and Susanne Kühnel (Jena, Germany) as co-workers during CT imaging sessions. We thank Benedikt Hallgrímsson (University of Calgary, AB, Canada) for providing time on the Scanco μCT35. T.K. expresses his gratitude to Tom Daniel (University of Washington, Seattle, WA, USA), who gave a helpful introduction into finite element analysis and provided access to the finite element software. Two anonymous reviewers provided valuable comments that helped to improve a previous version of this paper.