Predicting bite force in mammals: two-dimensional versus three-dimensional lever models

Measurements and estimates of whole-organism performance have provided useful information for studies of adaptation and phenotype–environment correlations. In vertebrates, bite force is an important performance trait that can be linked to whole-organism performance because it is relevant to several functions that may impact fitness (Anderson et al., 2008). Vertebrates use their cranio-mandibular apparatus not only to capture, subdue and process prey, but also in competition for mates (Korff and Wainwright, 2004; Lappin and Husak, 2005; Motta and Kotrschal, 1992; Reilly and Lauder, 1990). As a consequence, the ability to generate high bite forces can be crucial for determining the spectrum of available prey (Aguirre et al., 2003; Freeman and Lemen, 2007; Herrel et al., 2002b; Herrel et al., 2001; Huber et al., 2005), disputing and acquiring territories (Herrel et al., 1999; Lailvaux et al., 2004; Lappin and Husak, 2005; Vanhooydonck et al., 2005) and for males to win contests (Huyghe et al., 2005; Lailvaux et al., 2004).

In vivo bite force is well documented in reptiles (Erickson et al., 2003; Herrel and Holanova, 2008) but has been reported for only a few species of mammals (Aguirre et al., 2002; Dumont et al., 2009; Hylander et al., 1992; Santana and Dumont, 2009; Williams et al., 2009). Therefore, researchers have turned to models to estimate bite forces based on approximations of the structure of the skull and the physiology of the muscles that adduct the jaws. Bite force estimates are often used in comparative studies, where they are correlated with ecological tasks (Christiansen, 2007; Kiltie, 1982), stress distribution in the skull (Christiansen and Adolfssen, 2005; Thomason, 1991), and morphological variables such as bite point, gape, skull size or muscle mass (Dumont and Herrel, 2003; Herrel et al., 2008; Williams et al., 2009).

Thomason's method of predicting bite force is the model that is most commonly applied to mammals (Christiansen and Adolfssen, 2005; Christiansen and Wroe, 2007; Ellis et al., 2008; Thomason, 1991; Thomason et al., 1990; Wroe et al., 2005). The values of bite force predicted by this method are based on a static 2-D lever mechanics model in which muscle forces are determined from an assumed muscle stress and an estimated muscle area. These predicted values are consistently lower than in vivo measurements (Christiansen and Wroe, 2007; Ellis et al., 2008). It is not clear why this is the case, although one possibility is that the method by which muscle areas are determined is a source of error (Christiansen and Wroe, 2007).

Muscle force magnitudes are often calculated by multiplying muscle stress by muscle area. A common measure of muscle area is physiological cross-sectional area (PCSA) which can only be gathered from detailed dissections of fresh or preserved muscles. Not surprisingly, these data are available for only a handful of mammal species (Herrel et al., 2008; Perry, 2008; Santana et al., 2010; Taylor et al., 2009). Nevertheless, the widespread interest in bite force and its implications for organismal performance has lead to development of methods for estimating PCSA. One of the most common techniques used to estimate PCSA of jaw adductors was developed by Thomason (Thomason, 1991). Expanding on work by Kiltie (Kiltie, 1982; Kiltie, 1984), Thomason estimated PCSA from measurements of the infratemporal fossa in dry skulls as seen in 2-D images from posterodorsal and ventral views (Fig. 1). The area of the infratemporal fossa in the posterodorsal view was used to estimate PCSA of the temporalis muscle. In a ventral view of the skull, the area of the infratemporal fossa served as an estimate for PCSA of the combined masseter and medial pterygoid muscles.

Locations and directions of muscle forces relative to the TMJ axis affect the moment generated by each force. Kiltie (Kiltie, 1984) used linear distances collected from the sagittal view of the skull (the anteroposterior length of masseter scar) and mandible (the height of the coronoid process above the condyle) to approximate the placement of the muscle forces from which moments were calculated. Thomason's (Thomason, 1991) method places each force at the centroid of the area used to estimate PCSA. However, it is not clear that this is a reasonable approximation of the force centroid (i.e. the location at which a single force vector could be placed to create the same moment as the distributed force generated by a muscle). To incorporate more detailed 3-D anatomical information over which muscle forces are distributed, we developed a method for calculating bite forces using 3-D reconstructions of skulls and mandibles coupled with information about muscle attachment regions on the skull and mandible that were documented during specimen dissections. Previous methods have assumed that muscle forces act perpendicular to lines described in 2-D views of the skull (Kiltie, 1984; Thomason, 1991), but here we propose a method that utilizes the 3-D geometry of the skull and muscle attachment areas over which muscle forces are distributed. This 3-D method also aims each muscle force directly toward its respective insertion region on the mandible.

In this study, we test the hypothesis that bite force predictions generated from models that incorporate 3-D skull morphology and muscle anatomy from dissections provide more accurate bite force predictions than 2-D methods. Specifically, we predict that regressions of in vivo versus predicted bite forces based on our 3-D model will exhibit slopes closer to 1.0 and a higher r² value than those based on a 2-D model. We test this hypothesis and prediction using data from 24 species of New World leaf-nosed bats (family Phyllostomidae) and one species from a closely related family (Noctolionidae). We chose these species because of the availability of in vivo bite force data, muscle attachment regions identified during dissections and PCSA for the primary jaw adductors. Moreover, their wide range of skull morphologies provides a test of the broader utility of the methods for mammals with very different skull shapes.

In vivo bite force was measured for 24 phyllostomid and one closely related non-phyllostomid bat species (Noctilio albiventris). Bats were collected using mist nets at localities in Venezuela (2006, 2007), Panama (2007) and Mexico (Dumont et al., 2009). Only adult males and adult non-pregnant, non-lactating females were used for this study. Shortly after capture, we measured the bilateral molar bite force of each bat using a piezoelectric force transducer (Kistler, type 9203, range ±500N; Amherst, NY, USA; accuracy 0.01–0.1 N) attached to a handheld charge amplifier (Kistler, type 5995) (Herrel et al., 1999). The force transducer was mounted between two bite plates as described and illustrated in Herrel et al. (Herrel et al., 1999). The tips of the bite plates were covered with medical tape to protect the bats' teeth and to provide a non-skid surface. We adjusted the distance between the bite plates for each individual to accommodate a gape angle of about 30 deg. (Dumont and Herrel, 2003). We recorded at least five bite force measurements for each bat. Only unidirectional bite forces were measured. Shearing forces at the tooth were not measured with this device. Individuals for which we had complete bite force data were killed using cardiac compression, preserved in 70% ethanol, and kept as voucher specimens. All procedures were approved by the Institutional Animal Care and Use Committee at the University of Massachusetts, Amherst, USA (protocol # 26-10-06).

Several variables are required to calculate bite force with Eqn 1. These include muscle forces usually determined by multiplying muscle physiological cross-sectional area by muscle stress, the location of the bite point relative to the TMJ axis, and the locations of muscle forces relative to the TMJ axis. Although not used in bite force predictions based on the methods of Kiltie (Kiltie, 1982) and Thomason (Thomason, 1991), our 3-D model also incorporates the locations at which the muscles insert into the lower jaw. We use these muscle insertion sites to determine the directions of the muscle forces applied to the skull. This type of geometric information has been included in other models of the jaw, including those of reptiles (Cleuren et al., 1995; Herrel and Holanova, 2008), fish (Herrel et al., 2005) and bats (Herrel et al., 2008).

The second method of determining physiological cross-sectional area – what we call estimated PCSA – estimates PCSA using the 2-D techniques described by Thomason (Thomason, 1991), except that we used 2-D images of digitally reconstructed surface models in place of photographs of dry skulls. The skulls of the voucher specimens were cleaned using a dermestid beetle colony and scanned using a microCT-scanner (Skyscan 1172 Microfocus X-radiographic Scanner, Skyscan, Belgium) at Amherst College, Amherst, MA, USA. Details of scanning parameters are available from the authors upon request. The X-ray projection images produced by the microCT-scanner were converted using filtered back-projection into a volume consisting of a stack of X-ray attenuation cross-sections, or slices, using reconstruction software (NRecon v. 1.5.1.4, MicroPhotonics Inc., Allentown, PA, USA). For each species, the slices were imported into Mimics v. 13.0 (Materialise, Ann Arbor, MI, USA), where a thresholding tool was used to isolate the range of grayscale values representing the skull bones. Finally, a 3-D surface model consisting of contiguous triangles defining the shape of each skull was produced and saved as a stereolithograph (STL) file.

The estimated PCSA of the temporalis muscle was based on the area of the infratemporal fossa as seen in a posteriodorsal view of the skull (Fig. 1A). The plane of the posteriodorsal view was established using four points on the skull, the most lateral points of the left and right orbital processes and most posterior points on the left and right zygomatic arches. The left and right infratemporal fossae were outlined and their areas and area centroids (Fig. 1) calculated using ImageJ v1.42 (National Institutes of Health, Bethesda, MD, USA). The combined masseter–pterygoid estimated PCSA was estimated from the area of the infratemporal fossae as seen in a ventral view of the skull (Fig. 1). The ventral plane was defined by three points: the most distal tip of the premaxilla between the first pair of incisors, and the two occipital condyles. The regions occupied by the masseter and pterygoid muscles were outlined according to the method of Thomason (Thomason, 1991), and their areas and centroids (Fig. 1B) calculated using ImageJ v1.42.

In biomechanical models, muscle force is typically estimated by multiplying muscle area (measured or estimated PCSA) by muscle stress. Muscle stress has been shown in mammals to range between 147 and 500 kPa (Rohrle and Pullan, 2007; Sellers and Crompton, 2004; Thomason et al., 1990; Wilson and Bowers, 1975), and it is unlikely that there is a single value that is applicable to all species, or even to all muscles groups within an individual. One study has suggested that variation in muscle stress is associated with muscle function [i.e. flexor or extensor (Buchanan, 1995)] and we know that stress changes as muscle fibers (and sarcomeres) deviate from the optimum point on their length–tension curves (Anapol and Herring, 1989; Burkholder and Lieber, 2001; Gans and Devree, 1987; Meyer et al., 1998; Rassier et al., 1999; Rome and Lindstedt, 1997). Optimal fiber and sarcomere lengths were not collected for this study. Therefore we assumed a constant value of muscle stress of 250 kPa that is commonly used in the literature (Herrel et al., 2008; Nigg and Herzog, 1994). In our models we included the major muscles which act to close the jaw (temporalis, masseter, medial and lateral pterygoid), and modeled bilateral molar biting at a single gape angle of 30deg. Our models also assumed that all of the muscles were maximally activated.

We applied muscle forces to the skull in two different ways: either as point loads or distributed loads. When modeled as a point load (Point Load method), all of the muscle force is modeled as a single force vector acting perpendicular to the view from which PCSA is estimated and acting through its area centroid (Fig. 2). Scaling the magnitude of the muscle force to the value found using estimated PCSA exactly repeats Thomason's 2-D method (Thomason, 1991). In all of our models the bite point was estimated to be at the center of the first molar on the right side of the skull (and because of symmetry, the first molar on the left side, as well). The 3-D coordinates of the bite point were collected and used to calculate the vector from the TMJ axis to the bite point (i.e. r_j). We calculated bite force using Eqn 1. Additionally, by modeling bilateral biting we avoided the need to address different levels of muscle activation on working and balancing sides of the jaw (Hylander et al., 1992; Vinyard et al., 2006).

The second method of applying muscle forces to the skull was a 3-D Distributed Traction method that we developed to distribute forces over each muscle's respective attachment region on the skull as determined from dissections. We manipulated the 3-D surface models of the skulls using Geomagic v. 11 (Geomagic, Research Triangle Park, NC, USA) to open the jaw to a gape angle of 30deg. We then used digital images of the muscle dissections to define the 3-D regions of muscle origin on the skull and muscle insertion regions on the mandible for the temporalis, masseter, medial pterygoid, and lateral pterygoid (Fig. 3). Muscle origin and insertion regions were exported from Geomagic as 3-D surface models (STL files). Area centroids of the 3-D insertion regions were calculated using a Matlab (Mathworks, Natick, MA, USA) program (Area_Centroids_From_STL, available upon request). The inputs required for this program are the binary STL of the insertion regions in Cartesian coordinates, and the output is a spreadsheet file containing the Cartesian coordinates for the area centroid of each muscle insertion region.

An additional Matlab program (Boneload, available upon request) was written to distribute forces and calculate the moment contribution from each muscle group considering the 3-D geometry of the muscle attachment regions. Two files are necessary for running Boneload: a NASTRAN file containing nodal coordinates and connectivity for each of the muscle attachment areas, and a spreadsheet (in Excel *.xls format) containing individual muscle force magnitudes, insertion region centroid coordinates, and coordinates defining the TMJ axis. The NASTRAN file was created by importing STLs of muscle attachment areas into Strand7 (Strand7, Sydney, Australia), from which NASTRAN files were exported.

Assuming left–right symmetry and bilateral molar biting, we modeled only the right side of each skull. Muscle forces, calculated using muscle stress and PCSA, were distributed over the surface of each attachment region. We began by applying a uniform pressure on each STL triangle within the muscle attachment region on the skull. This pressure was equal to the total muscle force divided by the total area of the muscle attachment region. Next, we accounted for muscle stacking (Grosse et al., 2007) by adjusting the magnitude of the distributed forces within a muscle group. The force applied to each STL triangle within a muscle group was calculated as the uniform pressure multiplied by each individual STL triangle area. This force was scaled according to the ratio of the distance between the muscle insertion centroid and the furthest STL triangle of the muscle group to the distance between the muscle insertion centroid and the centroid of the element to which the force was applied. This method restricts the muscle forces on the skull that are furthest away from the insertion region (where there is little muscle stacking) to be lowest, and the forces that are nearest the insertion region to be highest. After accounting for muscle stacking, the magnitude of the total force applied to the attachment region was verified against the force calculated from muscle stress and PCSA. If necessary, force magnitudes were adjusted linearly across all the muscles until the total force applied matched the force calculated from muscle stress and PCSA. Finally, each muscle force was directed toward the area centroid of its respective insertion region on the mandible (Fig. 3).

In order to predict bite force, we summed moments about the TMJ axis due to muscle forces using Boneload. Although we assumed bilateral symmetry and modeled muscle forces and bite forces on the right side of the skull, the TMJ axis was defined as the line connecting the area centroids of left and right glenoid fossae. The bite force was calculated using Eqn 1 and assumed the bite force vector was perpendicular to the vector from the TMJ to the bite point. Shearing forces at the bite point were not considered because they could not be measured by the bite force meter. The bite force was then doubled and compared with bilateral molar bites collected from live animals in the field. We used ordinary least squares regressions of measured bite force (predictor variable) against calculated bite forces (response variable) to determine the slopes, intercepts, and coefficients of determination. We also compared measured PCSA against estimated PCSA of the temporalis muscle and masseter–pterygoid muscle combination using an ordinary least squares regression. Significance of the regression parameters was tested using t-tests. All statistical analyses were performed using SPSS for Windows (SPSS Inc., Chicago, IL, USA).

If all aspects of muscle function and skull geometry were accounted for correctly, there would be a perfect match between predicted and measured bite force. In statistical terms, the regression of measured bite force against predicted bite force should return a slope and an r² of 1.0 and an intercept of zero. Results from the regression of measured bite force against bite forces calculated using the Point Load method with estimated PCSA (Fig. 4), equivalent to the method developed by Thomason (Thomason, 1991), indicate a significant relationship between the measured and calculated bite forces (P<0.001; r²=0.613). The 95% confidence interval (CI) around the slope (β₁=0.858, 95% CI=0.563–1.15) encompasses 1.0 and the intercept does not differ significantly from zero (β₀=2.13, P=0.106). In contrast, while the regression of measured bite force on calculated bite forces using the Point Load method with measured PCSA (Fig. 3) is also significant (P=0.006), the 95% CI does not encompass 1.0 (β₁=0.263, 95% CI=0.085–0.442), the intercept differs significantly from zero (β₀=4.69, P=0.007) and the r² is low (0.288).

Of the two ways to generate PCSA values, the Point Load method with estimated PCSA provides reasonably accurate bite force predictions. However, it under-predicts measured PCSA for the temporalis muscle as indicated by a slope that is greater than 1 (β₁=1.72, Fig. 5). This method also over-predicts measured PCSA for the masseter–pterygoid muscle combination as indicated by slope that is less than 1 (β₁=0.637, Fig. 5). These results indicate that although the Point Load method with estimated PCSA predicts bite force well, it does so at the cost of misrepresenting the measured PCSA values.

The 3-D Distributed Traction method also demonstrates a strong association between measured and predicted bite forces. The regression of measured bite force against bite force predictions using the 3-D Distributed Traction method with estimated PCSA (Fig. 6) has an r² of 0.66 and an intercept that does not differ significantly from zero (β₀=2.07, P=0.087). However, the slope is greater than 1 (β₁=1.66, P<0.001, 95% CI=1.151–2.174). By contrast, the regression of measured bite forces against those predicted using the 3-D Distributed Traction method with measured PCSA, has a slope very close to unity (β₁=0.991, P<0.001, 95% CI=0.674–1.308). In addition, the intercept does not differ significantly from zero (β₀=2.33, P=0.057) and the r² value is relatively high (0.65). Overall, these analyses show that the 3-D Distributed Traction method with measured PCSA provides the closest match to the expected slope, intercept and r² values.

This study illustrates that the most accurate predictions of bite force are made with 3-D models that incorporate details of 3-D skull geometry, muscle attachment regions and dissection-based measures of muscle PCSA. Although 2-D methods provide a reasonable estimate of bite force, they inaccurately represent the anatomy of the muscles and fail to include the complex 3-D geometry of the skull and mandible in bite force estimates. If a researcher is interested only in estimating the magnitude of bite force, Thomason's 2-D method is an easily accessible approach that can provide reasonable estimates of bite force.

Although our 3-D model provided the best predictions of bite force, the method was not perfect. One assumption we made is that the subjects were in fact biting maximally during the collection of in vivo bite forces. We do not (nor may we ever) know this to be a fact, and behavioral variation in bite force could explain some of the variation in the data. Electro-physiological stimulus is one method of assuring that muscles contract maximally to produce a maximum measured bite force and has been done, for example, in domestic dogs (Ellis et al., 2009). Testing our 3-D Distributed Traction method against the data collected for domestic dogs would be one way to further investigate the accuracy of our 3-D model and validity of the assumptions we made. However, this method may not actually stimulate the entire muscle and would be extremely difficult to implement for small wild mammals in the field (Ellis et al., 2008).

Another potential source of variation, common to many studies that model jaw mechanics, is our understanding of muscle stress and its variation among species and muscle groups. We assumed a constant value of muscle stress in all of the adductor muscles and obtained reasonable regressions of measured bite forces against calculated bite forces. Since muscle stress modulates bite force linearly in these lever mechanics models, increasing or decreasing muscle stress will shift predicted values right or left on the x-axis of Figs 4 and 6. Regression intercepts will change, but slope and r² values, will remain unchanged. For this reason, the intercepts we report must be interpreted cautiously.

If, as seems likely, muscle stress varies among species or muscle groups within individuals, 3-D lever models that minimize geometric and anatomical approximations (e.g. skull morphology and PCSA), could be used to incorporate the effect of varying muscle force through force–length relationships of muscle fibers and sarcomeres, as has been done with jaws of fish (Van Wassenbergh et al., 2005; Van Wassenbergh et al., 2007). However, our results indicate that the assumption of a constant muscle stress in the 3-D Distributed Traction method with measured PCSA provides improved estimates of bite force over the 2-D methods with estimated PCSA. Allowing for muscle stress to vary among muscle groups according to force–length relationships could further improve the results. Herrel et al. (Herrel et al., 2008) recently demonstrated a significant correlation between in vivo bite forces, temporalis muscle mass and skull length in bats. We can use our 3-D model to further explore how these and other anatomical variables collectively modulate bite force. There are essentially only two variables that can be adjusted to change moments about an axis of rotation. They are the force vector, F_j, and the locations at which forces are applied relative to the axis of rotation, r_j (Eqn 1). The 3-D Distributed Traction method with measured PCSA can be used to investigate exactly how variables such as muscle mass and skull length collectively affect bite force. This advantage stems from replacing estimates of key variables with those directly measured from 3-D imaging and dissections (i.e. 3-D skull geometry, PCSA, muscle origin and insertion locations). Previous bite force models have included 3-D coordinates of origin and insertion regions (Cleuren et al., 1995; Herrel et al., 2002a; Herrel et al., 2008), thereby orienting muscle forces towards anatomically relevant locations. However, muscle force is modeled as a single force vector in those models. We present a novel method in which the muscle force is distributed over a 3-D muscle attachment area on the skull in addition to being oriented towards the insertions of the muscles on the lower jaw. With an accurate model of the biomechanics of bite force in hand, we can turn our efforts to investigating how key variables interact to produce bite force in living species, as well as to understanding their ecological specializations (Santana et al., 2010).

Fig. A1 illustrates the variables needed to predict bite forces when modeling either the mandible or skull as a third-class lever. Generally, 2-D methods of bite force predictions assume forces such as the masseter force, F_Mass, temporalis force, F_Temp, and bite force, F_Bite, act as point loads located at distances L_Mass, L_Temp and L_Bite respectively, away from the axis of rotation (i.e. TMJ axis).

This work was supported by a grant from the NSF (DBI 0743460 to E.R.D. and I.R.G.), a University of Massachusetts Natural History Collections David J. Klingener Endowment Scholarship and a Smithsonian Tropical Research Institute Predoctoral Fellowship to S.E.S. The authors would also like to thank Jaime Tanner and Tom Eiting for their discussions concerning this paper, Whitey Hagadorn for his technical assistance and allowing the use of the micro-CT scanner at Amherst College, and Dan Pulaski for his technical assistance with software.

Anatomy

 
• l_f

  muscle fiber length

 
• m

  muscle mass

 
• PCSA

  physiological cross-sectional area

 
• TMJ

  temporomandibular joint

 
• θ

  muscle pennation angle

 
• ρ

  muscle density

Mechanics

 
• F_j

  force vector

 
• F_Bite

  bite force vector

 
• F_Mass

  masseter force vector

 
• F_Temp

  temporalis force vector

 
• F_Bite

  magnitude of bite force vector

 
• F_Mass

  magnitude of masseter force vector

 
• F_Temp

  magnitude of temporalis force vector

 
• L_Bite

  distance from TMJ to bite force vector

 
• L_Mass

  distance from TMJ to masseter force vector

 
• L_Temp

  distance from TMJ to temporalis force vector

 
• j

  muscle index

 
• M_j

  moment about origin of r_j created by F_j

 
• M_j^TMJ

  magnitude of moment about the TMJ axis created by F_j

 
• n

  number of muscles

 
• r_j

  position vector locating origin of force vector

 
• TMJ

  TMJ axis vector

 
• ϕ

  angle between horizontal and temporalis force vector

 
• •

  vector dot product

 
• ⊗

  vector cross product

Statistics and regression

 
• CI

  confidence interval

 
• P

  significance value

 
• r²

  coefficient of determination

 
• β₀

  regression intercept

 
• β₁

  regression slope

Other

 
• STL

  stereolithograph