Local palaeoecological implications

In addition to a recently discovered procolophonid jaw [84], the Vellberg rhynchocephalian material represents the first small vertebrate remains from the source locality. It is generally agreed that Diphydontosaurus-like rhynchocephalians fed on small invertebrates [15, 91, 138, 147–149]. This is supported by their general body size, tooth shape and build of the lower jaw. The teeth are the same shape as tools that can puncture soft materials with relative ease, but are vulnerable to extreme torsion and bending [138, 147]. The slender jaws provide long out-levers for rapid closure and capture of small active prey but are not suited to withstanding substantial loading forces [149, 150]. The morphology of the Vellberg rhynchocephalian is consistent with these attributes, which is noteworthy as a predator of small invertebrates has not previously been described from this locality. By contrast, all of the taxa currently known from Vellberg (at least as adults) are suited to feeding on small vertebrates or fish. In turn, the Vellberg rhynchocephalian would itself have been prey for other animals in the community, such as immature individuals of Batrachotomus and other archosauromorphs. Like small vertebrates in modern communities (e.g. [151]), the lepidosaurs were probably an important link in the food chain between primary and tertiary consumers.

Global importance of the locality

The Middle Triassic record of small gracile vertebrates is poor. There are several rock units from around world that preserve terrestrial vertebrate remains: the Manda beds of Tanzania, Africa (e.g. [152–154]); the oldest part of the Santa Maria Formation (Fm) of Brazil (e.g. [155, 156]); the Chañares Fm of Argentina (e.g. [157]); the Moenkopi Fm of North America [158, 159]; the Yerrapalli Beds of India (e.g. [160, 161]); the upper part of the Beaufort Group of the Karoo Basin (e.g. [162, 163]); the Kelamayi Fm, Ermaying Fm, and Hongyanjing Fm of China (e.g. [164, 165]); the Donguz and Bukobay of Russia (e.g. [64, 166]); the Omingonde Fm of Namibia [167]; part of the Fremouw Fm of Antarctica [168]; the Zarzaïtine Series of Algeria [169]; the Areniscas y Lutitas del Figaro unit of Spain [170]; and a few units in the United Kingdom such as the Otter Sandstone (e.g. [142, 143]). However, fossils from the associated localities typically represent medium or large vertebrates such as trematosaurids, rhynchosaurs, cynodonts, and early archosaurs (e.g. [154, 160, 165]). Animals of small size (skull length <30 mm long) such as procolophonid reptiles are occasionally recovered but these are typically robust remains (e.g. [142, 155]). Therefore, as a new microvertebrate locality, Vellberg is expected to provide a more balanced picture of the Middle Triassic fauna and palaeoecological communities.

Vellberg may also shed light on the early fossil record of important tetrapod groups such as frogs, salamanders, caecilians, albanerpetontids, and choristoderes. All of these groups should have representatives in the Middle Triassic but currently none are known (e.g. [171–178]). Whether this absence of data is related to a failure to sample appropriate facies or a tendency for these animals to be small and gracile, or both, remains unclear [171]. Nevertheless, the material described here demonstrates that Vellberg has the potential to yield remains of other small tetrapods and to provide important information on a poorly known period of significant change in global ecosystems.

Divergence estimates and congruence with the fossil record

As an independent test of the internal ages of Squamata, we compared them to eight well described and dated fossils that could have been used as additional age constraints. All of them support our dating (Additional file 6), being as old or older than the mean of the estimate. However in three cases they would have truncated the younger bound of the credibility intervals by about 10 Mya.

Beside the manual control of eight alternative calibrations, we also ran the fossil-based cross validation analysis implemented in the penalized likelihood (PL) method of Near and Sanderson [131] on the 14 fossils used (Additional file 6). Simplified, the cross validation procedure sequentially removes one fossil at a time and estimates the node it constrains, to test whether a fossil causes a significant shift towards an older age of the node. Although this is not necessarily a problem with well described fossils, it may indicate a significant rate change close to that node that needs to be calibrated. Cross validation of our data set indicates that the most influential fossil is the calibration of crown-group Diapsida CNY (increased by 52 Mya, fraction score 0.24). This result is not unexpected as it is the node that constrains the root of the phylogeny (Additional file 6). The only other fossil that increases the age estimate significantly in the PL analysis is the fossil constraining the crown-group Anguimorpha CN12 (Figure 4), with an increase of ~8.6 Mya and a fraction score of 0.07.

Prior distribution of fossil-constrained ages

The setting of prior distributions for constrained ages is a non-trivial task. For the final dating analysis we used an exponential prior calibration density on divergence times. This approach means that the likelihood for the age of a node is highest at the age of the fossil, whilst the older possible ages have lower likelihood. Statistically the first ancestor of a lineage is not the oldest fossil known or recognisable clade member based on clear autapomorphies [179]. Therefore the use of the exponential prior is suboptimal, and the inferred ages are likely to be more prone to underestimation compared to other alternative priors. The widely used lognormal prior (see e.g. [39, 43]) implies that it is more likely for the real age to be older than the fossil. It can be argued that this prior would better represent the ghost lineage that must exist. However, in most cases there is no objective way of choosing the shape of the prior distribution, especially not in the case of organisms that are rarely preserved as fossils, and the analysis would potentially be highly biased toward the authors’ opinion on the fossil ghost range. Even if the lognormal prior could potentially approximate the true ages better, we chose the exponential prior because the minimum age of the fossil constraint is the only known date, the prior will be less biased toward the opinions of the researchers, and it represents a sound method from a philosophical viewpoint; our hypothesis is easily falsified if older fossils are found.

For comparison we also ran the analysis with less informative priors, such as a uniform prior with minimum age from the fossil and maximum age 50 million years older, which is probably too old for at least most of the constrained nodes. This approach tended to push all nodes towards older ages, resulting in median ages close to the given maximum age. It seems unlikely that all fossil records have a “ghost range” of almost 50 million years, so we do not recommend this prior distribution unless there is no other way of constraining the node.

The mean of the exponential prior was set consistently to 4.0 (in absolute ages a possible interval of 15–20 Mya with a small probability of estimated dates being older) to provide a plausible and hopefully useful dating of the lepidosaur and squamate divergence times. As a comparison we ran analyses with mean 5.0 and 6.0 (intervals of about 20 and 25 Mya, respectively). In general the resulting ages were within the range of ±5 My compared to the analysis with mean 4, but in some cases with larger credibility intervals (data not shown).

To check the influence of the highly informative priors on ages and monophyly, we also ran the same analysis but with sampling from priors only (by creating an empty alignment). The topology obtained was different, despite monophyly enforced on 16 nodes (beside the nodes constrained by fossil ages, mammals were set as sister to the reptiles (including Aves), and Squamata was set as monophyletic). The posterior density of ages on constrained nodes was similar for some of the nodes, whilst others had a very different distribution. This indicates that the results were not determined by the prior assumptions alone, but that the sequence data also influenced the results for all nodes.

Previous age estimates of lepidosaurs and squamates

Previous studies attempting to estimate the timing of the origin of Lepidosauria and Squamata have suggested a range of different ages (Table 1, Figure 5). The “outliers”, unusually recent or old dates, may be due to the method used, the number and selection of taxa, choice of molecular marker and of course the choice and number of fossil calibrations.

Figure 5 The phylogenetic relationships and fossil record of early lepidosaurs compared to molecular divergence estimates. Estimates for the origin of Lepidosauria based on previous molecular studies are listed on the right in blue with short arrows. Estimates for the origin of crown group Squamata are listed on the right in red with long arrows. Timescale based on Gradstein et al. [47]. Fossil records include those described, or referred to, in Butler et al. [105], Carroll [27], Clark and Hernandez [31], Colbert [30], Evans [8, 9, 26, 33, 91], Evans and Białynicka [34], Evans and Jones [5], Evans et al., [18], Fraser [22, 23, 136], Fraser and Benton [11], Heckert et al. [24], Nesbitt [180], Renesto [137], Reynoso [19, 150], Robinson [29], Sues and Hopson [13], Sues and Olsen [12], Whiteside [15], and others listed in Evans et al. [181] and Jones et al. [10]. Full size image

Gorr et al. [40] used a global clock approach to estimate divergence times within reptiles (including Aves). They concluded that there were large differences in evolutionary rates of reptilian hemoglobins between larger groups, causing an erroneous topology, so their age estimates should be viewed with caution. In a study on vertebrates, Kumar and Hedges [35] estimated gene-specific substitution rates, dated the separate gene trees, and then averaged over the trees to get one dated tree. As rates vary among lineages and therefore do not obey a global clock, they first excluded genes with extensive heterogeneity, and excluded the outliers before averaging over gene trees. Despite this, all nodes closer to the root showed large ghost ranges between estimates and first fossil record (e.g. Agnatha originating in the Precambrian), which is likely due to the method being unable to distinguish between extensive time or fast substitution rates. They conclude that the molecular ages are not overestimating the divergences, and that there are substantial gaps in the fossil record. Their estimate for Lepidosauria was 276±54.4 Mya.

Janke et al. [48] used mitochondrial genomes of a total 35 species to obtain rough estimates of divergence times for Squamata and turtles, assuming a constant evolutionary rate between 2 reference points: the Synapsida/Diapsida split (310 Mya) and the Crocodylidae/Aves split (254 Mya). Previously published genomes of a snake and side-necked turtle were excluded on the basis that their fast evolutionary rates complicate the phylogenetic analysis. The phylogeny and dating only contained two squamates (Iguana and the mole skink Eumeces). Neverthless, the origin of “Squamata” is stated to be 294 Ma and this date is argued to be consistent with the fossil record. However, this value actually represents the estimated divergence time between the lineage that includes the two squamates from one that includes turtles, crocodiles, and birds (Sphenodon was not used). Therefore this estimate more correctly represents the divergence time of Lepidosauromorpha rather than Squamata.

In general, studies using the Bayesian “multidivtime” method [182] give older age estimates than our study [37, 38, 49–51]. This is most likely an artifact of the method, which due to a strong autocorrelation assumption (the “minab” prior) tends to smooth ages towards the root of the tree to be consistent with the greatest tree depth (e.g. [36, 183]). This bias is also more prominent in studies with few taxa, and in the studies listed above the number of squamates range between 19 and 38. In practical terms this means that most of the employed minimum age constraints towards the leaves are uninformative in these “multidivtime” analyses. For the age of Lepidosauria we find no overlap between our deepest credibility interval (251.4 Mya) with the shallowest confidence intervals of Kumazawa [37] or Albert et al. [38], and in all the above mentioned studies using multidivtime crown-group Squamata is estimated to be Triassic or older.

Wiens et al. [41] used the semi-parametric penalized likelihood (PL) method of Sanderson [102] and 11 fossil constraints. For the PL analysis it is necessary to set a fixed age close to the root. The focus of the study was the internal divergences in Squamata, and they therefore chose to use the oldest known rhynchocephalian fossil to fix the most recent common ancestor of Squamata and Rhynchocephalia to 227 Mya. This approach may have caused an underestimation of the age of crown-group Squamata (178.7 Mya compared to 193.1 in this study). Hugall et al. [36] used RAG1 sequence and the PL method to study tetrapod diversification, including a total of 35 squamates. They used a maximum age of 450 Mya for the lungfish-tetrapod root and tested different sets of calibrations for internal nodes. All employed constraints were fixed, to avoid the method artefact of “model overfitting”, meaning that constraints closer to the leaves can lead to overestimation of deeper nodes. They estimated the split median ages between Sphenodon and Squamata to be 250–275 Mya, and that of crown-group Squamata to be 171–201 Mya, depending on the calibration scheme employed.

Hipsley et al. [42] used the same constraint for the lizard-tuatara split as Wiens et al. [41] but the former used a Bayesian probabilistic method as implemented in the TreeTime software [184]. To account for the uncertainty in fossil calibrations and the likelihood of the true age of a node being older than the first fossil record, the age constraint was set with a hard upper bound of 228 Mya and soft lower bound of 239.4 Mya. Their estimate for the Sphenodon-Squamata split was 238±10 Mya.

Pyron [39] proposed a method that can objectively test fossil placement and the likelihood of age estimates by comparisons between datasets of different studies. The empirical example is divergence analyses on RAG-1 DNA from 129 gnathostome taxa to compare the affect of two different sets of fossil calibrations. The sample included Sphenodon and 44 squamates. The uncorrelated lognormal method in BEAST was used, and a lognormal distribution was chosen for the prior distribution of ages from the fossil calibrations. Four fossil calibrations from Müller and Reisz [52] provided a mean estimate of 236 Mya (credibility interval 212–253) for Lepidosauria and 189 Mya (163–213) for Squamata whereas five fossil calibrations from Hugall et al. [36] provided a mean estimate of 265 Mya (240–290) for Lepidosauria and 208 Mya (179–234) for Squamata. The shallower estimates were preferred based on a comparison to the wider fossil record. These dates are similar to our own but have greater confidence intervals.

Mulcahy et al. [43] estimated divergence dates for squamates using 64 ingroup species and 25 nuclear loci (19,020 base pairs in total), comparing the results obtained from Penalized Likelihood (r8s) and the uncorrelated lognormal method in BEAST. The overlap between their study and the present one is substantial for terminal taxa. There are however some important differences in the fossil constraints such as the use of a younger rhynchocephalian fossil here (for a detailed comparison see Additional file 1). Mulcahy et al. [43] fixed the topology to the same maximum likelihood tree they used as input in the r8s analysis, to facilitate direct comparisons of ages between PL and BEAST, whilst we only constrained the calibrated nodes to be monophyletic. As opposed to our approach of using exponential age priors, Mulcahy et al. [43] applied lognormal distribution of ages for the 11 fossil constrained internal (Lepidosauromorpha) nodes. The oldest rhynchocephalian was set to 222.8 Mya, based on the Vinita specimen from the Ladinian–Carnian boundary. Note that this age was chosen using the timescale of Gradstein et al. [46] rather than the more recent Gradstein et al. [47]. The lognormal priors were set to have a mean and standard deviation of 1.0 – meaning a very narrow interval (about 3 Mya) with an arbitrary mean close to the minimum age of the fossil (e.g. for Lepidosauria 223.4 Mya, 222.9-225.9). Mulcahy et al. [43] conclude that the BEAST/lognormal clock analysis gives younger ages than the r8s/PL analysis. This is not surprising, considering that the internal priors have soft lower bounds but are strong enough to behave as if they have a hard bound, thereby constraining other internal nodes more than the minimum age constraints in the PL analysis, where the only lower bound is the fixed root. This is also likely to be the reason why the BEAST estimates seem more stable with narrower credibility intervals.

Origin time of Lepidosauria, crown-Squamata, and other major clades

The Vellberg jaw helps to bridge an important gap in the fossil record and establish that Lepidosauria (stem group Rhynchocephalia, and stem group Squamata) diverged at least 240 Mya (Figure 4). Discovery of lepidosaur remains in the Middle Triassic is consistent with previous predictions made by palaeontologists (e.g. [9, 26, 34]). It is also consistent with the Late Triassic rhynchocephalian fossil diversity [9, 23–25]. This new record from Vellberg supersedes previously used molecular dating calibration points of 223, 227, or 228 Mya for the lizard-tuatara split [37–43, 45].

For the divergence between Lepidosauromorpha and Archosauromorpha (bird-lizard split) the median of our estimate, 271 Mya (259–285), lies close to the boundary between the Lower and Middle Permian. This date is deeper than the oldest known fossils of either group (Figure 5): the earliest known archosauromorph is Protorosaurus from the Upper Permian (Wuchiapingian) of northeast England (UK) and the Kupferschiefer of Germany [104] and the earliest certain lepidosauromorph is Sophineta from the Lower Triassic (Olenekian) of Poland [5, 34]. However, the possibility that large gaps in the fossil record remain, particularly so for stem-lepidosaurs (>20 Mya), highlights the need to survey further fossil localities in the Middle and Late Permian for small vertebrates.

Importantly, our estimate strongly suggests that the origin of Lepidosauria postdates the Permian mass extinction event (252 Mya), which represents a significant period of environmental upheaval possibly linked to a runaway green house environment [65, 67, 68]. An Early-Middle Triassic origin and radiation of Lepidosauria would be associated with general changes from fairly uniform warm-arid environments towards ones experiencing humid-arid fluctuations and monsoon systems [59, 61, 62, 90]. Complex biodiversity was still in the process of reestablishment after the Permian end mass extinction event [65, 67, 68]. Vegetation in the Middle Triassic was dominated by gymnosperms such as cycads, ginkos and conifers [59, 62]. Coeveal macrofaunal changes include the diversification of early archosaurs such as the sail-backed poposaurs and appearance of the first dinosauriformes (e.g. [105, 152, 180, 185]). The subsequent “Carnian Pluvial Event” (CPE) of the Late Triassic [61] is thought to represent a global increase in rainfall and further shifts towards more humid climates (Figure 4).

Our results suggests that the origin of crown-group Squamata lies in the Early Jurassic,190 Mya (175–212). We cannot exclude the possibility that crown-squamates appeared before the late Triassic extinctions but our median estimate post-dates them. Our estimate lies soon after changes in general vegetation that indicate changes towards warmer climates and greater continental aridity [54, 62]. This may be part of a general shift towards more regionalised climates and environments, at least in the northern hemisphere. The earliest secure lizard fossils currently referred to crown-Squamata are Middle Jurassic in age and therefore do not conflict with our estimate [8, 17, 18].

The Cretaceous origin of most major crown-groups suggests the radiation of Squamata occurred after and alongside continental fragmentation (Figure 4, Table 3, [58, 186]. Therefore the widespread distribution of many modern lizard groups today (e.g. [1]) probably requires a number of post Jurassic dispersal events to have occurred. Evidence that transoceanic dispersal of squamates is possible does certainly exist (e.g. [187–190]) and the distances between continental fragments in the second half of the Mesozoic were much less than they are today [191].

Our estimates for the origin of most modern groups coincide with a general improvement of the squamate fossil record [5–8] and contraction of rhynchocephalian distributions to southern continents [10, 20, 181, 192]. This shift in lepidosaur communities may be related to expansion/contraction of preferred environments [6, 7, 21] or displacement by active competition [20, 192, 193], but distinguishing between the two hypotheses remains problematic [5, 25]. The Early Jurassic to Early Cretaceous diversification of crown-group squamates is concurrent with that of several modern lineages of beetles [54, 56]. However, rather than reflecting a predator–prey relationship it may be symptomatic of the general development of more modern ground cover and microhabitats.

The divergence estimates for both crown-group Iguania, 136 Mya (117–152), and total group Iguania 153 Mya (148–161) post date the fossil taxon Bharatagama from India originally referred to Iguania [18]. It is possible that Bharatagama represents an early stem crown-group squamate with a jaw morphology convergent with modern acrodont iguanians, or that it belongs to another clade.

Our estimated origin time for Gekkota, 76 Mya (52–101), is younger than that of some previous studies but there is some overlap between credibility intervals (e.g. [36, 43, 194]). There are also two early-mid Cretaceous fossils that could potentially challenge our crown-group age of Gekkota: Cretaceogekko burmae preserved in amber from Myanmar (>97.5 Mya) [195] and Hoburogekko suchanovi from Mongolia (125–99.6 Mya) [196, 197]. Both fossils likely belong to the gekkotan lineage but their precise relationship with extant geckos is unclear [197]. Cretaceogekko was described as crown-group gekkotan based on it’s advanced adhesive toe pads, but it has recently been inferred that specialized toe pad morphology has evolved (and been lost) several times across the gekkotan phylogeny [198]. Hence it is not possible from morphological characters alone to determine crown-group affinity. The redescription by Daza et al. [197] of Hoburogekko concludes that the combination of jaw and skull characters is likely to belong to a stem-group gekkotan, and that a phylogenetically conservative placement of these Cretaceous fossils is recommended.