Jaw biomechanics and the evolution of biting performance in theropod dinosaurs

Manabu Sakamoto


Despite the great diversity in theropod craniomandibular morphology, the presence and distribution of biting function types across Theropoda has rarely been assessed. A novel method of biomechanical profiling using mechanical advantage computed for each biting position along the entirety of the tooth row was applied to 41 extinct theropod taxa. Multivariate ordination on the polynomial coefficients of the profiles reveals the distribution of theropod biting performance in function space. In particular, coelophysoids are found to occupy a unique region of function space, while tetanurans have a wide but continuous function space distribution. Further, the underlying phylogenetic structure and evolution of biting performance were investigated using phylogenetic comparative methods. There is a strong phylogenetic signal in theropod biomechanical profiles, indicating that evolution of biting performance does not depart from Brownian motion evolution. Reconstructions of ancestral function space occupation conform to this pattern, but phylogenetically unexpected major shifts in function space occupation can be observed at the origins of some clades. However, uncertainties surround ancestor estimates in some of these internal nodes, so inferences on the nature of these evolutionary changes must be viewed with caution.

1. Introduction

Theropod dinosaurs display a high degree of craniomandibular morphological diversity, which can lead one to infer equally diverse feeding mechanisms, and thus feeding preferences, strategies, behaviours and ecology (Bakker 1986; Paul 1988; Henderson 2000; Barrett 2005; Barrett & Rayfield 2006). Biomechanical models have been used as a quantitative means to investigate the mechanics of feeding in theropods (Molnar 2000; Henderson 2002; Therrien et al. 2005), and such models have reached a new level of complexity with the advent of high-resolution computer modelling methods (Rayfield 2007). Such comparative biomechanical studies rely on the assumption, even though rarely explicitly stated, that the selected taxa truly represent different function types (such as ‘typical’ or ‘specialized’), but the presence and distribution of function types across Theropoda have rarely been assessed. Thus, we are left with an unanswered fundamental question: how do biting performances differ among theropod taxa? Here, a simple but novel method of ‘biomechanical profiling’ is devised to capture overall biting performance along the entirety of the tooth row. Further, multivariate ordination is employed to visualize ‘function space’ (Anderson 2009) occupation.

Another important question that is seldom asked is: what is the evolutionary history or the underlying phylogenetic structure of biting performance in theropods (Westneat 1995, 2004; Sakamoto et al. 2010)? Without a proper phylogenetic framework, biomechanical studies lack evolutionary context. There are numerous methods available for phylogenetically informed comparisons (phylogenetic comparative methods, PCMs; Harvey & Pagel 1991), but it is essential first to test for phylogenetic signal because PCMs may not produce reliable results if the data are not correlated with phylogeny (Laurin 2004). Testing for phylogenetic signal in itself is also of interest, as this allows further interpretation about the evolution of biting performance in theropods.

2. Material and methods

(a) Biomechanical profiling

Biting performance is quantified using a specifically defined and well-established biomechanical metric, mechanical advantage (MA; Westneat 1994), the ratio between the moment arm of the muscle moment (effort) and the moment arm of the biting moment (load). MAs were computed for each biting position along the entirety of the tooth row (figure 1) in 41 theropod taxa and the outgroup Plateosaurus, using a simple two-dimensional lever model. Potentially valuable three-dimensional information can be lost by this simplification, but the benefits of using explicitly defined simple models outweigh the drawbacks; for example, their explicitness and statistical rigour. Muscle reconstructions are based on anatomical observations of modern avian and crocodilian jaw adductor musculature (M. Sakamoto 2005–2008, personal observation; but also Holliday & Witmer 2007; Holliday 2009; electronic supplementary material, figure S1).

Figure 1.

Biomechanical profiling by application of a simple biomechanical model and standardization of tooth row length. The skull of Sinraptor (modified from Currie & Zhao 1993) is shown with attachments and lines of action of three major muscle groups; M. adductor mandibulae externus group (MAME), M. adductor mandibulae posterior group (MAMP), and M. pterygoideus group (MPT). The distances from the jaw joint (J) to these lines of action are the effort arms (mE, mP and mPt for each muscle group, respectively). Load arms mBi are taken as the distances between J and bite points along the entirety of the tooth row. Biting positions are represented as percentages of the tooth row, with 0% and 100% being the posterior- and anterior-most positions, respectively.

Since theropods have vastly differing tooth counts, direct comparisons are not possible without standardization of biting positions. Following their computation from absolute distances, MAs were plotted against corresponding biting positions scaled as percentage of tooth row length (figure 1 and electronic supplementary material, figure S4). Then a polynomial function of either the second order (y = β0 + β1x1 + β2x22), third order (y = β0 + β1x1 + β2x22 + β3x33) or fourth order (y = β0 + β1x1 + β2x22 + β3x33 + β4x44) was fitted (electronic supplementary material, figure S4) based on the fit of the nth-order polynomial determined by Akaike information criteria weights in the paleoTS library (Hunt 2008) in R (R Core Development Team 2009). The three common coefficients of the polynomial functions (β0, β1 and β2) were subjected to principal components analysis (PCA) using PAST (Hammer et al. 2001) to visualize function space occupation by theropod clades.

(b) Phylogenetic comparative methods

The phylogeny of Lloyd et al. (2008) formed the basis of all comparative analyses (electronic supplementary material, figures S5 and S6).

(i) Testing for phylogenetic signals

In order to test for phylogenetic signals in the biomechanical variables, the phylogenetic eigenvector regression (PVR; Diniz-Filho et al. 1998) was employed (Sakamoto et al. 2010). PVR is based on a multiple linear regression in which phylogeny, the predictor variable, is represented as eigenvectors (or principal coordinates axes) of a phylogenetic distance matrix. A multivariate multiple regression was performed in R with the polynomial coefficients as the response variable matrix and the first 34 eigenvectors (explaining 95% of phylogenetic variance) as the predictor variables.

Another test for detecting phylogenetic signal devised by Blomberg et al. (2003) was conducted using the picante library (Kembel et al. 2009) in R. This method uses phylogenetically independent contrasts (PIC; Felsenstein 1985) and compares the variances of the PIC computed from a given variable on a particular tree topology with those computed from permutations of the variable across the same tree. If the variances in PIC for the real data are lower than those from the permutations, then there is a significant phylogenetic signal (Blomberg et al. 2003). Blomberg et al.'s (2003) K statistic was also computed; K < 1 would indicate that closely related taxa have values that are less similar than expected under Brownian motion evolution (departure from Brownian motion, such as adaptive evolution), while K > 1 would indicate that closely related taxa have values more similar than expected, thus a strong phylogenetic signal (Blomberg et al. 2003).

(ii) Tracing the evolution of function space occupation

To trace the evolution of function space occupation, nodal values (ancestors) were estimated for the three polynomial coefficients separately using the maximum-likelihood (ML) method of ancestor character estimation (Schluter et al. 1997) available in the APE library (Paradis et al. 2004) in R. PCA was conducted on ancestor estimates and terminal values to visualize the evolution of function space occupation across phylogeny.

3. Results

(a) Biomechanical profiling

Biomechanical profiling reveals a continuum of various MA profiles along the vertical axis (figure 2a). The basal taxa Plateosaurus and Herrerasaurus are situated at the higher end of the vertical axis. Many taxa share a similar profile, with a gradual parabolic decline in MA values from the 0–100% tooth row positions, the only major difference being their vertical positions. However, some taxa exhibit profiles with steeper slopes and parabolic curvatures, while some show the opposite with shallower slopes and parabolic curvatures.

Figure 2.

Biomechanical profile plots and function space occupation. (a) Best-fit polynomial lines for MA against per cent tooth row in 41 theropod taxa and Plateosaurus represent biomechanical profiles of biting efficiency along the entirety of the tooth rows. Vertical position relates to the magnitude of MA, slope relates to the degree of change in MA along the tooth row and curvature relates to the degree of change in slope along the tooth row. (b) A plot of the first two principal components (PC) scores of the polynomial coefficients visualizes function space occupation. PC1 and PC2 account for 89.9% and 9.9% of the variance, respectively (both axes are to scale). The relationships of the polynomial coefficients with the PC axes are also shown. Tetanuran taxa form a continuous distribution, with a major axis corresponding to the spectrum of biting efficiency. Labels: 1, Plateosaurus; 2, Herrerasaurus; 3, Dilophosaurus; 4, Syntarsus; 5, Coelophysis; 6, Ceratosaurus; 7, Carnotaurus; 8, Majungasaurus; 9, Dubreuillosaurus; 10, Eustreptospondylus; 11, baryonichine; 12, Spinosaurus; 13, Monolophosaurus; 14, Sinraptor; 15, Yangchuanosaurus; 16, Allosaurus; 17, Acrocanthosaurus; 18, Carcharodontosaurus; 19, Gorgosaurus; 20, Daspletosaurus; 21, Tarbosaurus; 22, Tyrannosaurus; 23, Dilong; 24, Guanlong; 25, Compsognathus; 26, Ornitholestes; 27, Proceratosaurus; 28, Garudimimus; 29, Gallimimus; 30, Ornithomimus; 31, Struthiomimus; 32, Erlikosaurus; 33, Citipati; 34, Kamyn Khondt oviraptorosaur; 35, Ingenia; 36, Khaan; 37, Archaeopteryx; 38, Sinornithosaurus; 39, Dromaeosaurus, 40, Velociraptor; 41, Bambiraptor; 42, Deinonychus.

PCA of the polynomial coefficients distinguishes these differences (figure 2b); there is a clear distinction between the Tetanurae (but also Carnotaurus and Majungasaurus) and the basal taxa Plateosaurus, Herrerasaurus and coelophysoids (but also Eustreptospondylus, Ceratosaurus and Carcharodontosaurus). Taxa with steeper slopes and curves score higher along PC1 while taxa with shallower slopes and curves score lower. Separation along PC2 arises mostly from the contrast between vertical positioning and parabolic curvature. Tetanuran taxa are distributed along a major axis inclined at roughly 25° from the PC1 axis (figure 2b) with the ‘high-efficiency’ function types towards the top half and ‘weak/fast’ function types towards the bottom half of the distribution.

(b) Phylogenetic comparative methods

(i) Testing for phylogenetic signals

A multivariate PVR test shows a significant correlation between the polynomial coefficients and the first 34 phylogenetic eigenvectors (tables 1 and 2). Blomberg et al.'s (2003) test also detects significant phylogenetic signals in all three coefficients and results in K > 1 (table 3).

View this table:
Table 1.

Overall significance of the multivariate phylogenetic eigenvector regression. Approximate F-values, degrees of freedom and p-values are shown for each test statistic. ***Significant at 0.001.

View this table:
Table 2.

Univariate significances of the regressions of each response variable. Multiple R2, adjusted R2, F- and p-values are presented for each response variable. **Significant at 0.01; ***significant at 0.001.

View this table:
Table 3.

Strength and significance of phylogenetic signal in biomechanical coefficients using Blomberg et al.'s (2003) method. K-statistic indicating strength of phylogenetic signal, observed variance of the PIC (varobs), mean variance of PIC from the permutation (varrnd), Z-score of observed versus random variance of PIC and p-value of observed versus random variance of PIC (Kembel et al. 2009) are shown for each variable (transformed coefficients). ***Significant at 0.001.

(ii) Tracing the evolution of function space occupation

PCA of the polynomial coefficients and their ML ancestor node estimates reveals that the common ancestors of Saurischia and Theropoda (figure 3, nodes 43 and 44, respectively) occupied a similar position in function space. While the common ancestor of Neoceratosauria (node 49) retains close to the ancestral function space, there is a shift towards the lower spectrum of the PC1 axis at the common ancestor of Tetanurae (node 51), with subsequent clades within Tetanurae all originating further along this axis. While Allosauroidea (node 57) originate at the high-efficiency end of the tetanuran distribution, Coelurosauria (node 61) and Tyrannosauroidea (node 66) originate close to the weak/fast biter end (figure 3). Both Ornithomimosauria + Maniraptora clade (node 70) and Maniraptora (node 74) originate at the weak/fast end of the tetanuran distribution, but the common ancestor of Therizinosauroidea and Oviraptorosauria (node 75) plots out towards the high-efficiency end.

Figure 3.

Evolution of function space occupation. The evolution of function space occupation can be traced through ML ancestor estimation and PCA on the terminal and nodal values. PC1 and PC2 account for 90.6% and 8.9% of the total variance, respectively (axes are not to scale; PC2 is exaggerated for graphical purposes only and should not be read literally). Internal nodes are displayed as grey points and nodes at major clades are indicated by arrows. Internal and terminal branches are shown in solid and dotted lines, respectively. Labels: 43, Saurischia; 44, Theropoda; 48, Neotheropoda; 49, Neoceratosauria; 51, Tetanurae; 56, Avetheropoda; 57, Allosauroidea; 61, Coelurosauria; 62, Tyrannosauroidea; 63, Tyrannosauridae; 70, Maniraptora + Ornithomimosauria; 74, Maniraptora; 75, Oviraptorosauria + Therizinosauroidea; terminal labels as in figure 2.

4. Discussion

(a) Biomechanical profiling and function space occupation

Biomechanical profiling using MA along the entirety of the tooth row allows for investigation of functional disparity patterns across Theropoda (figure 2). Biomechanical profiles not only quantify the magnitudes of MA, but also subtle mechanical characteristics of tooth row morphology in relation to the effort arms. The vertical position of the profile (β0) relates to the proximity of the tooth row to the jaw muscles, while the slope (β1) relates to the length of the tooth row. The parabolic curvature (β2) relates to a combination of these two characteristics; steeper curves are associated with long tooth rows that are in close proximity to jaw muscles.

Function space occupation by theropod taxa shows interesting disparity patterns. At the high-efficiency end of the tetanuran function space distribution, tyrannosaurs, allosaurs and ceratosaurs plot out close to each other despite their obviously different cranial morphology. Allosaurs have long been considered to be weak biters (Bakker 1998; Rayfield et al. 2001), but strikingly, their MA profiles are here found to be similar to those of abelisaurid ceratosaurs. Although MA does not take into account the magnitude of muscle forces, previous reconstructions of muscles however have also yielded similar bite forces for Allosaurus and Carnotaurus (Rayfield et al. 2001; Mazzetta et al. 2004, 2009), further supporting functional similarities in these two disparate clades. An extreme case of high-efficiency biting and functional convergence with ceratosaurs can be observed in Carcharodontosaurus, which plots beyond the range of distribution of tetanuran theropods and converges with Ceratosaurus towards the ancestral function space.

Derived tyrannosaur profiles, on the other hand, are distinguished from allosaurs and ceratosaurs in having shallower slopes. In particular, Tyrannosaurus exhibits higher consistency in MA along the entirety of its tooth row (figure 2), enabling it to maintain high bite force more rostrally than in other theropods. Oviraptorosaurs have lower β0 values, so were less powerful biters than other high-efficiency biters, but their extremely high β1 values make them the most efficient at the 100 per cent tooth row position.

Weak/fast biters exhibit rostrally displaced tooth rows (high mB) along with posteriorly positioned/oriented muscles (low mM). Their consistently low MA values throughout their tooth rows indicate they were less efficient at delivering forceful crushing bites, but more effective in fast snapping bites. Weak/fast biters are typically small- to medium-sized theropods, with the exception of spinosaurs; fast biting is suitable for piscivory. Dromaeosaurs are here found to be low-efficiency biters, but this does not preclude the possibility that they were ‘big-game hunters’ (Paul 1988, p. 364) because if they did indeed use their sickle claws, then ‘the jaws were secondary weapons’ (Paul 1988, p. 364). Further, weak-biting Varanus komodoensis (Moreno et al. 2008) is capable of inflicting deep wounds through saw-motion biting (Auffenberg 1981), so it is not difficult to imagine dromaeosaurs similarly delivering saw-motion bites if kicks were not enough to kill their prey.

(b) The evolution of biting performance in theropods

The strong phylogenetic signal in polynomial coefficients (tables 13) can be observed in phylogenetic reconstructions of ancestral function space occupation (figure 3); with a few exceptions, most theropods plot out in function space relatively consistent with phylogeny or following a phylogenetic trajectory. It is also noteworthy that, while changes close to the terminals are relatively small (except for coelophysoids), there are major shifts at the bases of clades, in particular Neoceratosauria (node 49), Tetanurae (node 51), Allosauroidea (node 57), Tyrannosauridae (node 63) and the Oviraptorosauria + Therizinosauroidea clade (node 75). This implies either (i) that these large changes arise from methodological artefacts of ML ancestor estimation, (ii) that taxonomic (and temporal) coverage close to the origins of clades (basal members or sister taxa) are not sufficient, or (iii) that evolution of biting performance across Theropoda actually occurred more rapidly at the origins of clades. The former two implications are most certainly applicable for the basal nodes (especially node 48) in which ancestor estimates are quite likely affected by computational artefacts (figure 3). Similarly, tyrannosauroid sampling suffers from a long temporal gap between the basal and derived taxa, so a major shift in function space occupation at the base of Tyrannosauridae is not at all surprising considering branch duration. On the other hand, there is good support for the suggestion that the origin of the Oviraptorosauria + Therizinosauroidea clade was associated with a genuine shift in function space occupation; under Brownian motion evolution, this clade would be both phylogenetically and temporally (electronic supplementary material, figure S6) expected to originate in function space proximal to parent and sister nodes but also to basal coelurosaur taxa (figure 3). The same may be true also for Allosauroidea considering its relationship in function space with the common ancestor of Coelurosauria, but the tentative nature of ancestor estimates for Tetanurae and thus Avetheropoda warrants caution in this interpretation.

An important observation is that the coelophysoids Coelophysis and Syntarsus do not represent a primitive condition for Theropoda with respect to biting performance, but instead are unique functional specialists; their basal phylogenetic positions are not reflected in their biomechanical profiles.

(c) Phylogenetic structure of variance in biting performance

The presence of a significant and strong phylogenetic signal (tables 13) in theropod biomechanical profiles indicates that the evolution of biting performance mirrors phylogeny. Phylogenetic distance is the primary source of difference in biomechanical profiles, and unique selection pressures that may cause differential rates of change in individual branches are rare (non-departure from Brownian motion). This is a striking revelation since it counters the intuition that functional characters would be expected to be under unique selective pressures. However, the underlying mechanism of this phylogenetic signal, whether constraint or inertia (Blomberg et al. 2003; McKitrick 1993), may never be known.

(d) Issues with biomechanical modelling

In reality, musculoskeletal systems are complex, both geometrically and physiologically. Any form of biomechanical modelling necessitates a certain level of approximation. Different models allow for different degrees of sophistication: for instance, a dynamic versus a static model, two-dimensional versus three-dimensional or inclusion versus exclusion of various parameters. Dynamic models can model in vivo motions with reasonable accuracy (Moazen et al. 2008), and can also estimate static biomechanical values such as bite force with relative accuracy (Moazen et al. 2008). However, static models have also been demonstrated to estimate bite force accurately (Sinclair & Alexander 1987; Cleuren et al. 1995; Herrel et al. 1998). Differences in two- and three-dimensional models are yet to be investigated in detail, but a recent study would indicate that the two models result in comparable bite force estimates (McHenry et al. 2007). Inclusion of additional model parameters can certainly add sophistication but at the same time increase potential error. This is primarily because of the fundamental uncertainties associated with the estimation of these variables, many of which remain to be validated. Recruitment patterns of muscle activity are difficult to determine in vivo (Busbey 1989) and are approximated using electromyographical data (Cleuren et al. 1995), which are not available in fossil taxa. Nevertheless, inclusion of such information in models for extant taxa did not result in bite force estimates that were significantly different from those estimated by simple models assuming simultaneous muscular action (Cleuren et al. 1995). The effects of muscle pennation can be substantial in estimating physiological cross-sectional areas of muscles (Biewener & Full 1992), but is often complicated and difficult to quantify even in living forms (Busbey 1989; Hieronymus 2006; M. Sakamoto 2005–2008, personal observation). Again, there is no means of estimating muscle pennation in fossil forms. On the other hand, it is possible to estimate muscle size from some osteological correlate, but such estimates are simply a matter of scaling (i.e. reflecting size of the bony element). More importantly, size can span several orders of magnitude (compare Compsognathus and Tyrannosaurus), dwarfing any size-independent biomechanical variation; the model would effectively be capturing size and not function. Several authors have sought to eliminate this size-associated variance to isolate the biting function by either post hoc scaling (Wroe et al. 2005; Christiansen & Wroe 2007; Sakamoto et al. 2010) or a priori scaling (Dumont et al. 2009; Slater & van Valkenburgh 2009).

MA by definition quantifies a specific biomechanical performance and has been used extensively in previous studies (Busbey 1989; Westneat 1994, 1995, 2004; Anderson 2009; Sakamoto et al. 2010). The benefits of using a simple two-dimensional biomechanical metric are in its explicitness, limited uncertainties and statistical rigour. By taking a simple but explicit metric, it is possible to isolate important biomechanical properties from various forms of noise and address broad-scale questions such as functional disparity, evolution of biomechanical performance and the role of phylogeny.


The author thanks David Norman, Paul Jeffrey, Paul Barrett, Zhonghe Zhou, Fucheng Zhang, Makoto Manabe, Jean-Yves Sire, Hajime Taru and John Hutchinson for access to collections and specimens. Thanks also to James Rohlf, Fred Bookstein, Norman MacLeod, Joseph Felsenstein, Michel Laurin and Øyvind Hammer for discussions on multivariate and phylogenetic analyses. Many thanks go to Mark Bell for help with R, Graeme Lloyd for help with R (and especially for R codes: http://www.graemetlloyd.com/methdpf.html) and comments, and Mike Benton for proof reading this manuscript, but also for continued guidance and support. Thanks also to the two anonymous reviewers who helped refine this manuscript. This research was partly funded by the BBSRC.

  • Received April 14, 2010.
  • Accepted May 14, 2010.


View Abstract