Mathematical Model of the Human Jaw System Simulating Static Biting and Movements After Unloading

G.E.C. Slager1, E. Otten1, T.M.G.J. van Eijden2, and J. D. van Willigen1

1 Department of Medical Physiology, University of Groningen, 9712 KZ Groningen, The Netherlands; and 2 Department of Functional Anatomy, Academic Center for Dentistry Amsterdam, 1105 AZ Amsterdam, The Netherlands

    ABSTRACT
Abstract
Introduction
Methods
Results
Discussion
References

Slager, G.E.C., E. Otten, T.M.G.J. van Eijden, and J. D. van Willigen. Mathematical model of the human jaw system simulating static biting and movements after unloading. J. Neurophysiol. 78: 3222-3233, 1997. When the resistance to a forceful isometric bite is suddenly removed in unloading experiments, the bite force drops to zero and the mandible reaches a constant velocity. This occurs at an initial bite force of 100 N after ~12 ms when the incisors have moved 4.5 mm. Reflex activity is far too slow to limit the velocity at impact. To explore the influence of other factors (cocontraction, force-length properties, and force-velocity properties of the muscles) on the velocity at impact, a numerical forward dynamic model of the jaw system is formulated. Unloading experiments in different experimental conditions were simulated with the model. Most parameter values of the model are based on physiological data, both from literature and a data basis from a human cadaver study. Other parameter values were found by optimally fitting the model results to data from the unloading experiments. The model analysis shows that the limitation of the jaw velocity mainly may be due to the force-velocity properties of the jaw-closing muscles. Force-length properties of the jaw muscles hardly contribute to the impact velocity. The compliance of tendinous sheets in the jaw muscles is unfavorable for the reduction in impact velocity, whereas cocontraction of jaw-opening and -closing muscles helps to limit impact velocity. The force-velocity propertiesof the muscles provide a quick mechanism for dealing withunexpected closing movements and so avoid damage to the dental elements.

    INTRODUCTION
Abstract
Introduction
Methods
Results
Discussion
References

The human oral system is so powerful that it can easily damage its dental elements. Biting through hard and brittle food especially can result in high impact velocities of the dental elements. Any mechanism that reduces the bite force as soon as a movement commences may prevent damage to the teeth because the impact velocity will be lower.

In unloading experiments (in which the resistance to a forceful static bite is suddenly withdrawn experimentally), it is shown that the bite force does indeed decrease at high rate as soon as the mouth starts closing (Miles and Wilkinson 1982; Van Willigen et al. 1997). After the bite force vanishes after ~12 ms (Nagashima et al. 1997), the velocity no longer grows, because by virtue of Newton's law the velocity of a solid body is constant when no forces are exerted on the body. The aim of the present study is to uncover the factors that may cause this quick decrease in bite force and the limitation in the velocity of the mandible.

Reflex events can be excluded because these events occur too late (Hannam et al. 1968; Lamarre and Lund 1975; Miles and Wilkinson 1982; Van Willigen et al. 1997; Yoshida and Inoue 1995) to have any significant influence on control of the movement of the mandible.

Cocontraction of the jaw-opening and -closing muscles might be a factor that can cause a quick drop in force. Miles and Wilkinson (1982) suggest that the resistance to elongation of the activated cocontracting digastric muscles (resulting in an increase in jaw-opening force) is possibly responsible for the quick decrease in bite force and the limited mandibular velocity. They suggest that the resistance is perhaps due to distortion of cross bridges between myofilaments (Rack and Westbury 1974). However, Van Willigen et al. (1997) show that various levels of cocontraction have hardly any influence on the decrease in bite force during mouth closure.

Also the force-length properties of muscles can play a role in the drop in bite force. Slager et al. (1995) show in a model study (with linearized muscle models, not based on morphometric data) that the quick decrease in bite force can be attributed partly to the force-length characteristics of the jaw-closing and -opening muscles.

Furthermore, the force-velocity properties of the jaw muscles can have profound effects on the dynamic bite force (Otten 1991; Slager et al. 1995): the jaw-closing muscles can lose a fair amount of their force when they shorten, and the opening muscles can gain force when they are stretched, resulting in vanishing of the bite force.

However, from unloading experiments only, it is hard to judge how cocontraction, force-length, and force-velocity properties contribute to the decline in bite force during jaw movement (especially because the length and the velocity of the muscle fibers is unknown due to compliance of tendinous sheets). Numerical models can serve as a tool to find their relative contribution.

With this in mind, we have formulated a mathematical model that can simulate unloading experiments. By removing the force-length or force-velocity properties of the jaw muscles, the compliance of the tendinous sheets, or by changing the level of cocontraction of the jaw muscles, their contribution to the dynamic bite force and thus to the magnitude of the impact velocity can be studied.

    METHODS
Abstract
Introduction
Methods
Results
Discussion
References

We have formulated a mathematical forward dynamic model of the jaw system that can simulate jaw unloading experiments. For that we needed a model of the unloading device (device model) with which experiments were done; muscle models of the jaw-closing muscles (hereafter muscle model 1) and jaw-opening muscles (muscle model 2), based on morphological, physiological, and biomechanical properties; and to tune the model, we used results (force profiles, positions, and velocities of the mandible) from unloading experiments by Nagashima et al. (1997) (hereafter, the experiments).

The experimental set-up can be summarized as follows. An "unloading" device (Fig. 4A) was used comprising two parallel aluminum bars of which the upper was attached to two vertical plates mounted on a base plate. The lower bar hinged around an axis. Bite forces were exerted between the upper and lower incisors and cuspids on both bars of the device. For recording bite forces, strain gauges were attached on the lower bar. The initial resistance to closing was achieved by an empowered solenoid. The solenoid could be switched off---triggered by the output of the lower bar strain gauges---at a voltage equivalent to 100, 80, 60, and 40 N (hereafter, initial bite force). When the solenoid was switched off, the lower bar dropped at the back, and its bitten end was lifted up. The displacement of this bitten end was measured. The distance of travel of the lower bar could be varied by means of an adjusting screw. To buffer the shock of collision, this screw was covered with a rubber cap. To vary the initial mouth opening, the position of the upper bar could be adjusted by means of slotted holes and bolts.


View larger version (34K):
[in this window]
[in a new window]
 
FIG. 4. A: unloading device: 2 parallel aluminum bars (with pairs of strain gauges) are attached to 2 metal plates mounted on a base plate; the lower bar is fixed to an axis. Initial resistance to closing is achieved by a solenoid. Solenoid could be switched off at voltages equivalent to 100, 80, 60, and 40 N. B: an example of a model simulation of the unloading device at an initial loading force of 100 N and a distance of travel of the lower bar of 3.7 mm. It shows the counter force produced by the unloading device (Fresistance, thick line). Before unloading, Fresistance is equal to the initial loading force. After unloading, Fresistance rapidly declines, its profile being the sum of Foffset, Fmagnet, Fvibration, and Fend.

Five subjects were asked to bite without visual guidance through a resistance of 100, 80, 60, or 40 N with care. This was done at four initial mouth openings (24.5, 27.5, 30.5, and 33.5 mm) and three distances of travel of the lower bar (1.1, 2.4, and 3.7 mm), giving 48 different experimental conditions. The mouth opening (MO) of the subjects was defined as the interincisor distance added to the vertical overlap of the dentition (overbite). All experiments were repeated five times per subject. We had access to the raw data of these experiments.

Model of the jaw system

The model of the jaw system calculates forces, velocities, and positions of the mandible as a function of the initial bite force, mouth opening, and distance of travel of the mandible after the moment of unloading.

Figure 1 illustrates a diagram of the model of the jaw system. As can be seen, the model is built up from three compartments.


View larger version (33K):
[in this window]
[in a new window]
 
FIG. 1. Diagram of the model of the jaw system used for simulating the unloading experiments of Nagashima et al. (1997). Forces involved are coded as: 1, Fclosers; 2, Fresistance; 3, Fopeners. Unloading starts with a quick reduction of Fresistance (A). Result is an acceleration of the mandible (with massm) and the lower bar (with masslb) of the unloading device (see Fig. 4A), accompanied by a concentric contraction of the jaw-closing muscles and an eccentric contraction of the jaw-opening muscles (B).

MUSCLE MODEL 1. All the jaw-closing muscles (the masseter, temporalis, and medial pterygoid muscles), together were modeled as one single muscle (muscle 1), attached to the fixed skull and to the mandible. The mandible is attached to the lower bar of the unloading device by means of a stiff spring, simulating the suspension of the dental elements.

MUSCLE MODEL 2. The jaw-opening muscles (the digastric, mylohyoid, geniohyoid, and lateral pterygoid muscles) also were modeled as a single muscle (muscle 2), attached to the fixed hyoid bone and to the mandible.

DEVICE MODEL. The resistance of the unloading device was modeled by a suitable set of mathematical formulae, adequately describing the measured properties of the device (see Device model).

Apart from the intramuscular movements, the model of the jaw system has only one degree of freedom; i.e., the movement of the lower bar and mandible along the line of action of both muscles. [This is in contrast with the work of Laboissière et al. (1996), who have produced a model of the jaw system with 7 muscles with separate jaw and hyoid bone movements. Such a sophisticated model is useful in the context of the study of multimuscle control systems, but would be out of place in our more limited scope.]

The line of action of our model has a direction that is perpendicular to the occlusal plane and passes through the canines. (In reality the jaw muscles have a distributed attachment to the lower jaw, which has a variable center of rotation. To keep the model simple and yet realistic, we calculated the ratios between the force produced by the jaw muscles and the resultant force at the lower canines. The resultant force was used, because this way we could compare directly the external measured forces in dynamic equilibrium with the model output.)

The mentioned ratios are the kinematic transmissions from the muscles to the teeth. The actual muscles then could be described by two model muscles, which have properties that compensate for the simplified geometric arrangement of the model.

The forces produced by muscles 1 and 2 were called Fclosers and Fopeners. Because the net muscle force produced by opening and closing muscles together is consumed partly by acceleration of the mandible, the force measured on the lower bar (Foutput) is less than the total net muscle force. We therefore defined Foutput = Fclosers - Fopeners - massm · accm in which massm is the mass of the mandible and accm is the acceleration of the mandible. The resistance of the unloading apparatus as calculated by the device model was called: Fresistance and the initial bite force was called Fstart.

At the start of the simulation, the modeled muscles contract isometrically, exerting a bite force equal to Fstart, so that Fclosers - Fopeners = Fstart, and Fstart = Fresistance = Foutput. After time t = tunloading, Fresistance drops quickly, resulting in an acceleration of the lower bar (with masslb) and the mandible (with massm) accompanied by a concentric contraction of muscle 1 and an eccentric contraction of muscle 2 (Fig. 1B).

The model of the jaw system simulates a short static phase before unloading and the dynamic phase after unloading until the preset distance of travel is reached. In the simulations, the lower bar of the device and the mandible were described as separate solid bodies with masses. When we simulated the behavior of the oral system in free motion, we removed the mass and resistance of the device from the model.

Figure 2 offers a flow diagram of the simulations. The initial conditions of the model simulations were Fstart and Fopeners (see APPENDIX: tuned parameters), the initial mouth opening, and the travel distance of the lower bar.


View larger version (59K):
[in this window]
[in a new window]
 
FIG. 2. Flow diagram of the model. Input of the model are the initial conditions: the initial bite force (Fstart) and force of muscle 2 (Fopeners), the initial mouth opening, and the distance of travel. Output of the model are the bite force (Foutput), the velocity, and the position of the lower bar.

We calculated Foutput and the velocity and position of the lower bar as follows: Fclosers was calculated by Fclosers = Fstart + Fopeners. We derived Fresistance from the device model (see Device model). Because Foutput = Fclosers - Fopeners - massm · accm, we could define the driving force on the lower bar as Fdrive = Foutput - Fresistance. The acceleration of the solid body representing the lower bar could be calculated from this Fdrive and its mass. From this, by numerical integration, we were able to calculate the velocity and position of the lower bar and mandible.

Numerical integration with time steps of 0.1 ms allowed accurate calculation of the changing time-, position-, and velocity-dependent forces (Fclosers, Fopeners, and Fresistance). The recruitment of muscles 1 and 2 was calculated in the first time step from Fstart and Fopeners by inverting Otten's (1987a) formula for the activation dynamics of muscle fibers, producing muscle force from recruitment. The recruitment of both modeled muscles was kept constant throughout the whole simulation.

Muscle model

Muscle models 1 and 2 were founded on the muscle model written by Otten (1987a); a model of a tendinous sheet was added; Fig. 3 gives a diagram.


View larger version (29K):
[in this window]
[in a new window]
 
FIG. 3. Diagram of the muscle model. Muscle model contains a model of lumped muscle fibers, a model of a tendinous sheet, and a small internal solid body. Active force (Factive) is calculated from the recruitment, fiber length (curve A), and fiber contraction velocity (curve B). Passive force (Fpassive) is calculated from the strain of the passive structures of the muscle fibers (curve C). Force exerted by the muscle fibers (Ffiber) is equal toFactive + Fpassive. Force exerted by the tendinous sheet (Ftendon) is dependent on its strain (curve D).

Each muscle model contains a model of lumped muscle fibers, a model of a tendinous sheet, and a small internal solid body (with a mass) representing the effective moving muscle tissue.

The force exerted by the contractile part of the muscle model (Ffiber) is the resultant of active forces (Factive) [calculated from the recruitment, fiber length (curve A) and fiber velocity (curve B)] and passive forces (Fpassive) [depending on the strain of the muscle belly (curve C)]. The force exerted by the tendinous part of the muscle model (Ftendon) depends on the strain of the tendinous sheet (curve D).

Note that the length and the velocity of the contractile part of the muscle model depends on the position and velocity of the mandible and the length and rate of length change of the tendinous sheet.

Because the contractile part and the tendinous part of the model are positioned in series (Fig. 3), the force exerted by the muscle model equals Ffiber and equals Ftendon before unloading. After unloading, however, Ffiber - Ftendon becomes nonzero so that the internal solid body is accelerated. To avoid intramuscular oscillations, we added damping properties to the tendinous part of the model, in line with Hannam and Langenbach (1995).

Muscle models 1 and 2 differ in the values of their morphometric and muscle parameters (see Parameters used in the muscle model).

Device model

The purpose of the device was simply to quickly let the resistance to the bite decline in a reproducible way: no attempt was made to perform any servo tracking of the resulting movement.

The resistance offered by the device (Fresistance) after unloading appeared to be a complex ensemble of force components. After recording the resistance of the device in a large collection of tests, a search procedure was started to find the parameters of a set of mathematical formulas describing these force components (see Device parameters, APPENDIX).

Four force components can be discerned: Foffset (the resultant of dry friction between the lower bar and apparatus and an imbalance of the lower bar), Fmagnet (being the decaying remnants of the magnetic field), Fvibration (a time varying force originating in a damped vibration of the lower bar), and Fend (formed by the 1st contact of the lower bar with its rubber rest) apart from the force generated by accelerating the lower bar of the device. Figure 4B gives an example of a model simulation of Fresistance (thick curve) and its force components (initial loading 100 N; travel distance 3.7 mm).

Parameters used in the models

The parameters used in the muscle models can be subdivided into three sets: parameter values from literature (Table A1), morphometric parameter values taken from a cadaver study by one of us (van Eijden) (Tables A2 and A3), and parameters of which only rather wide ranges were found in literature. The values of the last set of parameters could be found by tuning of the model (Table A4A). For the value of the passive damping of the tendons, we have chosen a specific value (Table A4B). The parameters are described in the APPENDIX.

 
View this table:
[in this window] [in a new window]
 
TABLE A1. Parameter values from literature used in both muscle models

 
View this table:
[in this window] [in a new window]
 
TABLE A2. Two-dimensional morphometric data

 
View this table:
[in this window] [in a new window]
 
TABLE A3. Parameters based on morphometric data of a cadaver study

 
View this table:
[in this window] [in a new window]
 
TABLE A4. Parameters of the muscle models

Sensitivity analysis

To develop some notion on the relative meaning of model parameters for the force output, we performed a sensitivity analysis by testing the effect of an increase of 1% of the value of each parameter separately on the percent change in force output of the model.

Search procedure

After setting the morphometric parameters and the literature-determined parameters in the model, a search procedure was started to optimize parameter values by obtaining a minimal least squares fit of the model results on the experimental data (see APPENDIX).

    RESULTS
Abstract
Introduction
Methods
Results
Discussion
References

Comparing model output with experimental results

Figure 5 shows the fit of the model simulations (thick lines) with the experimental data (thin lines, SD in gray, n = 25) of eight experimental conditions. The panels depict forces, velocities, and positions of the lower bar at initial forces of 40 and 100 N, mouth openings of 24.5 and 33.5 mm, and travel distances of 1.1 and 3.7 mm. A comparison is made from the time of unloading until the preset distance of travel was met by the model (white time blocks).


View larger version (54K):
[in this window]
[in a new window]
 
FIG. 5. Results of 8 groups of experiments of Nagashima et al. (1997) and their model simulations. Depicted are forces, velocities, and positions of the lower bar (thin lines: experimental results with standard deviation in gray bands; thick lines: model results) at an initial bite force of 40 and 100 N, a mouth opening (MO) of 24.5 and 33.5 mm, and a distance of travel (TD) of 1.1 and 3.7 mm. Model results are only depicted in the white area's in each panel in the figure.

As can be seen, the model output is in good agreement with the experimental results, which also is illustrated in Fig. 6. The latter figure gives impact velocities as calculated by the model plotted against the measured average impact velocities of the experiments. As can be seen the data lie close to the line of equality.


View larger version (22K):
[in this window]
[in a new window]
 
FIG. 6. Impact velocities calculated by the model plotted against the averaged impact velocities as measured in the 48 experimental conditions of Nagashima et al. (1997).

Factors contributing to the reduction in bite force after unloading

Figure 7 gives information on the internal mechanics of muscle 1 and 2 in one simulation of an unloading experiment of 100 N (mouth opening 33.5 mm, travel distance 3.7 mm). It shows the time course of forces, lengths, and velocities of various components of muscles 1 and 2. The figure serves as an example. In the following text, numerical references will be made to other conditions that have not been included in Fig. 7.


View larger version (22K):
[in this window]
[in a new window]
 
FIG. 7. Data on the internal mechanics of muscles 1 and 2 in a simulation of an unloading experiment of 100 N (initial mouth opening = 33.5 mm, TD = 3.7 mm). Time course of various forces (i.e., the contribution of Factive, Fpassive, Fclosers, and Fopeners to Foutput), the lengths and the velocities of the components of muscles 1 and 2 are shown. A: note that the quick reduction in Foutput is caused mainly by the sharp decrease of Factive of muscle 1, representing the closing muscles. B and C: note that after unloading, the changes in length and the velocity of the muscle fibers is much smaller than that of the muscle.

Figure 7A illustrates for the jaw-closing (top) and jaw-opening muscles (bottom) the contribution of Factive, Fpassive, Fclosers, and Fopeners to Foutput. From Fig. 7A, top, it can be seen that before unloading Fclosers mainly consists of active muscle force (82%). [At a smaller initial bite force (40 N) at the same mouth opening, the active force of muscle 1 takes up a much smaller part of Fclosers (52%).] The passive force of muscle 2 (bottom) is almost zero, so that muscle 2 delivers active force only.

After unloading, the bite force (Foutput) drops from 98 to 9 N within 13 ms. Looking at the contribution of the different components, the active force of muscle 1 decreases rapidly from 84 to 4 N, whereas its passive force decreases from 19 to 15 N. At the same moment, the active force of muscle 2 increases from 5 to 12 N. Furthermore, 2 N is produced by deceleration of the lower jaw. Apparently, the quick reduction in Foutput is caused mainly by the sharp decrease of the active force of muscle 1 (closing muscles), whereas the changes in passive force of muscles 1 and 2 hardly contribute to the effect.

Figure 7, B and C, illustrates the relative length and velocity of muscles 1 and 2 and of their muscle fibers and tendinous sheets. The length of the structures is set to zero at the start of unloading to make comparison easy.

It can be seen that after unloading the changes in length and velocity of the muscle fibers are much smaller than that of the muscle. [At a change in length of 3.7 mm of muscle 1, the length change of the muscle fibers is only 1.8 mm (48%)]. This is due to the in-series arrangement of the tendinous sheet and the fibers and to differences in their mechanical properties. After unloading, the tendinous sheet and the fibers suddenly shorten with the tendon initially taking up most of the slack length (<= 52%). At the end of the movement, the length of the tendinous sheet increases again, whereas the muscle fibers still shorten. A similar effect, though more short lasting, can be seen in muscle 2.

Figure 8 shows the force-length (Fig. 8A) and force-velocity properties (Fig. 8B) of muscles 1 and 2. The ranges covered during a simulated experiment also are shown (conditions: 100 N initial force, 33.5 mm mouth opening and 3.7 mm travel distance). Black symbols indicate the start of the experiment, whereas white symbols illustrate its end. Normalized length is defined as the length of the muscle divided by its optimal length, which is the length at which the muscle produces its highest isometric active force.


View larger version (19K):
[in this window]
[in a new window]
 
FIG. 8. Influence of the active and passive force-length properties (A) and force-velocity properties (B) of the fibers of muscles 1 and 2 on the bite force in a simulated unloading experiment (initial bite force, 100 N; initial mouth opening, 33.5 mm; TD, 3.7 mm). A: it can be seen that, because of their force-length properties, muscles 1 and 2 showed a small increase in force during the dynamic phase after unloading. B: because of the force-velocity characteristics of the muscle fibers, Factive of the concentric contracting fibers of muscle 1 dramatically decreases whereas Factive of the eccentric contracting muscle fibers of muscle 2 increases marginally in its absolute value.

The figure illustrates that after unloading the fibers of muscle 1 shorten from a normalized length of 1.22-1.19, whereas those of muscle 2 lengthen from 0.71 to 0.78. (At a mouth opening of 24.5 mm and the same travel distance these figures are 1.10 to 1.08 for muscle 1, and 0.91-0.98 for muscle 2). The relative change in fiber length of muscle 2 is larger (7%) as compared with that of muscle 1 (2-3%). This is due to the fact that muscle 2 has a shorter length than muscle 1 and that there is less tendon creep of muscle 2 due to the muscles low activation.

It appeared that in all 48 experimental conditions simulated---due to their active force-length properties---muscle 1 as well as muscle 2 showed a small increase in active force when the force velocity effects are disregarded (the isometric active force). At an initial force of 100 N and a travel distance of 3.7 mm, the increase in isometric active force of muscle 1 is between 1.2 and 4.1 N and between 0.1 and 1.9 N of muscle 2. These values are too small to contribute significantly to the large reduction in bite force after unloading.

Note that the passive force-length properties of muscle 2 have no influence on the model results because muscle length was well below the rising part of the force-length curve.

From the force-velocity relationships (Fig. 8B), it can be seen that after unloading, Factive of the concentric contracting fibers of muscle 1 dramatically drops from 84 to 4 N at a mouth opening of 33.5 mm (and at a mouth opening of 24.5 mm from 97 to 9 N). The active force of the eccentric contracting muscle fibers of muscle 2 increases from 5 to 9 N at a mouth opening of 33.5 mm (and at a mouth opening of 24.5 mm from 5 to 8 N). Clearly, the force-velocity properties of muscle 1 (the closers) are mainly responsible for the steep decline in bite force after unloading in the simulations.

Factors contributing to the magnitude of the impact velocity after unloading

Figure 9 shows the impact velocities of the mandible in simulated unloading experiments of 40, 100, and 200 N. The initial mouth opening is set at 33.5 mm and the travel distances are 2.0, 4.0, 6.0, and 8.0 mm. In the simulations, the initial force of muscle 2 was set to be 5.1 N. 


View larger version (47K):
[in this window]
[in a new window]
 
FIG. 9. Impact velocities as predicted by the model for an unconstrained human mandible (no unloading device present), at an initial force of 40, 100, and 200 N, an initial mouth opening of 33.5 mm, and TD of 2.0, 4.0, 6.0, and 8.0 mm. Contribution of the various biomechanical factors to the impact velocity was studied by leaving them out of the model.

To simulate "natural" biting through hard and brittle food, the influence of the unloading apparatus was switched off (Fresistance = 0 N after tunloading and masslb = 0 kg). The histogram in Fig. 9A shows the results. As can be seen, the impact velocity does not increase at distances of travel >4-6 mm (a dynamic equilibrium is reached between the opening and closing forces at that point).

The contribution of the active force-length characteristics of the muscle fibers to the magnitude of the impact velocity was studied by leaving them out of the model (Fig. 9B). Without force-length properties, the impact velocity increases by 31% at a travel distance of 8.0 mm and Fstart = 40 N, by 13% at Fstart =100 N, and by 7% at Fstart = 200 N.

By leaving the force-velocity characteristics of the muscle fibers out of the model, the impact velocity increases tremendously (Fig. 9C); in the 40 N condition at a mouth opening of 33.5 mm and over a travel distance of 8.0 mm, an increase of 162% is calculated, and in the 100 and 200 N condition, both have an increase of 194%. Differences between the 40, 100, and 200 N simulations mainly are due to the passive muscle components which are not velocity dependent. In the 40 N experiments, the contribution of these passive components is relatively large.

The contribution of the compliance of the tendinous sheet of the jaw-closing muscles to the impact velocity was explored by reducing it to zero so that the sheet could not take up length (Fig. 9D). This gave rise to a reduction of impact velocity of 37% (40 N simulations), 63% (100 N simulations), and 70% (200 N simulations) at a distance of travel of 8.0 mm. Apparently, the compliance of a tendinous sheet in the jaw-closing muscles favors the impact velocity.

Finally, the influence of the cocontracting jaw muscles on the impact velocity was explored by increasing Fopeners from 5.1 N to 20 N (Fig. 9E). This increase in cocontraction caused the impact velocity to become 0.03 m/s after a distance of travel of 2.2 mm in the 40 N simulations. In the 100 N and 200 N simulations, the impact velocity levelled out at 0.06 and 0.16 m/s above a travel distance of 4.4 and 7.6 mm, respectively.

Sensitivity analysis

The sensitivity analysis shows a dominant influence of the parameters of the passive force-length curve of the muscle fibers in both muscles, the parameter determining the stiffness of the tendinous sheet of both muscles, the optimal length of the fibers and tendinous sheet of muscle 1, the fiber type composition of muscle 1, the initial force exerted by the cocontracting jaw-opening muscles, and the mass of the lower bar. An increase in parameter value of 1% of any of the above dominant parameters gives a change in the model force output of 0.1-0.9%. The model is robust to changes in the other parameters of the jaw model, giving a change in the model force output of 0.0-0.07% at a 1% increase.

    DISCUSSION
Abstract
Introduction
Methods
Results
Discussion
References

In this study, we showed that the model fits the experimental results offered by Nagashima et al. (1997) well. Not only do the impact velocities match, but more importantly, the time-varying output forces that result in these impact velocities match also in great detail. Figure 6 illustrates that the model predicts somewhat lower impact velocities than measured. This is due to noise in the position signal of the experiments; maximal velocities were calculated from this noisy position signal, resulting in slightly higher values of the impact velocities than would be available after curve fitting.

Limitations of the model

In modeling, one is forced to choose a compromise between realism and elegance. When a model is complex and data are far and few between, too many unknown values of parameters make the model easy to match any data set, whereas the underlying processes may not be in agreement with reality. At the opposite end, an elegant model may not be able to simulate the actual background of some observed phenomenon. In the present study, we have looked for a compromise. Therefore, we started with a model that was obviously too simple and added elements until the model was able to simulate the experiments. The good agreement between model and experiment is primarily due to the facts that sufficient subprocesses were modeled (such as intramuscular dynamics) as well as sufficient muscle physiological properties and that the model contains seven parameters (Table 4A) that could be tuned within sufficiently wide physiological boundaries. To keep the number of free parameters low, we looked for data in literature.

The composition of the model is formed by two dynamic skeletal muscle models connected in an antagonistic way to a representation of the lower jaw, which connects to a model of the experimental device. We founded the muscle models on the muscle model of Otten (1987a).

Basically these muscle models are a Hill-type model (Hill 1938), which has force-velocity properties that are not history dependent. This choice may give some unrealistic effects in force production, which may have been covered up by parameter tuning.

To come to anatomic realism, morphometric data from a cadaver study were used in the muscle models. Different action lines of the muscle fiber bundles of a single jaw muscle were simplified to one action line and the diversity in compartmentalization and variation in recruitment of motor units within a single muscle was ignored.

The geometric transformation of the coordinates of the sites of attachment of the opening muscles resulted in a rather small value of the length of muscle 2. This caused a relative large change in muscle fiber length of the openers at small movements of the jaw. However, the human hyoid bone moves considerably during normal oral behavior (Pancherz et al. 1986; Thexton et al. 1981). Because in the present model approach we did not wish to simulate hyoid movements, we kept the position of the hyoid bone constant. This invokes an overestimation of the muscle fiber velocities of the opening muscles. The hyoid bone had a different position from that found in the morphometric data set (see APPENDIX, Morphometric parameters) to avoid opening muscle insufficiency at large mouth openings.

Recruitment of both muscles was kept constant during the whole simulation, which is based on the observation that before unloading the bite force was relatively constant in the experiments and that after the unloading, reflex events evoked in the jaw muscles are too late to have an effect on the bite force during the dynamic phase (Miles and Wilkinson 1982; Van Willigen et al. 1997). This is in contrast with studies on perturbations of the limb, in which reflex events can have an influence on the control of the limb movement because the duration of the movement is much longer (Angel et al. 1965; Dufossé et al. 1985; Soechting and Lacquaniti 1988).

The optimization procedure found a value of 5.1 N for Fopeners. Comparable low levels of cocontraction are experimentally found for the anterior digastric muscles (Miles et al. 1982; Van Willigen et al. 1997) and for the inferior head of the lateral pterygoid muscle (Yoshida and Inoue 1995). Information on the activation of other opening muscles is unavailable.

Results of the model study

The model simulations suggest that the impact velocity of the mandible is limited by the force-velocity properties of the jaw-closing muscles. Without these force-velocity properties the impact velocity would increase by a factor of 1.94-2.7 and 3.9 m/s at an initial force of 100 and 200 N, a mouth opening of 33.5 mm and a travel distance of 8.0 mm; the kinetic energy at impact would be 8.6 times larger in this situation. We therefore conclude that the force-velocity properties of the jaw-closing muscles are a key factor in preventing damage to the teeth after sudden unloading of the jaw.

The active and passive force-length properties of the muscles hardly play any role in limiting the jaw impact velocity. The simulations suggest that small differences in impact velocities at different mouth openings (at equal initial forces and travel distances) are to be attributed to the passive force-length properties of the jaw-closing muscles. The relative contribution of these passive forces to the total force is greatest at large mouth openings and small initial forces. These passive forces do not decrease at increasing jaw velocity in contrast with active forces.

The model suggests that cocontraction of the jaw-opening muscles---if sufficiently high---helps in reducing the impact velocity. However, this factor may be of less importance than the simulations suggest due to the wrong assumption that the position of the hyoid bone is fixed during mouth closure. Thexton et al. (1981) and Pancherz et al. (1986) describe a movement of the hyoid bone in an upward-forward direction during mouth closure, which limits the velocity of the fibers of the opening muscles.

Model simulations with an uncompliant tendinous sheet of the jaw-closing muscle predict a reduction of impact velocity of 37-70% over travel distances of 8.0 mm. Apparently this compliance is unfavorable over these distances. In the simulations, the tendinous sheet takes up most of the muscle length change after unloading so that the velocity of the muscle fibers remains low, reducing the loss in active muscle force due to the force-velocity properties of the fibers. One could ponder on the functional meaning of tendinous sheet compliance in human jaw-closing muscles. Although it is unfavorable for the impact velocity, compliance may be an unavoidable design factor when packing many short fibers in a limited space because long tendinous sheets usually are encountered in this kind of muscle architecture.

The model achieves a good fit without using the effects of reflex events. This is in contrast with studies on perturbations of the limb, in which reflex events can have an influence on the control of the limb movement because the duration of the movement is much longer (Angel et al. 1965; Dufossé et al. 1985; Soechting and Lacquaniti 1988).

Looking at the results in a broader context, our study suggests that, in particular, the force-velocity properties of the active agonist muscles offer a powerful mechanism for handling effects of sudden perturbations. This mechanism is quick and independent of neural feedback. We submit that the physiological properties of muscles and tendinous sheets should have a more prominent place in theories on motor control.

    ACKNOWLEDGEMENTS

  The authors thank Drs. P.J.W. Jüch and M. L. Broekhuijsen for helpful comments regarding the text of this manuscript. The authors are willing to submit the computer source of the model described to anyone who applies for it.

    APPENDIX: PARAMETERS, USED IN THE MUSCLE MODEL

Literature-fixed parameters

Table A1 gives the literature-fixed parameters, used in both muscle models.

In concurrence with Otten (1987a), we characterized the force-velocity relationship of the muscle fibers by a constant k (determining the curvature of the relationship) and Vmax (the maximal shortening velocity of the muscle fiber). Both parameters are muscle fiber-type dependent and because the jaw-closing muscles differ in fiber composition from the jaw-opening muscles (Eriksson et al. 1981, 1982; Eriksson and Thornell 1983), we used a fiber-type composition-dependent value of k. This value was found by linearly interpolating between the values k = 0.17 for slow fibers and k = 0.25 for fast fibers (see Table A1). The interpolation was guided by the fiber type composition (see Tuned parameters). The same holds for Vmax, in which Vmax = 7.1 fiber lengths/s for slow fibers and Vmax = 18.3 fiber lengths/s for fast fibers. Both parameters are based on Close (1964). Vmax is the maximal shortening velocity at full activation. The shortening velocity attainable depends on recruitment (Julian and Moss 1981) and was calculated from Vmax and recruitment using a relationship from Otten (1987a).


View larger version (15K):
[in this window]
[in a new window]
 
FIG. A1. Calculated active force-length curve for human sarcomeres. Dashed line, curve produced by means of the sarcomere model of Otten (1987b) with data from Walker and Schrodt (1974) as input. Solid line, curve produced by a formula suggested by Otten (1987a) for the active force-length relationship. Curve B was fitted to curve A with coinciding maxima.

The active force-length relationship of the muscle fibers was described using a formula with three parameters, a, b, and s, which determine the roundness, skewness, and width of this curve (Otten 1987a).

Walker and Schrodt (1974) published values for the length of actin and myosin filaments of human muscle fibers. These values were fed into the sarcomere model of Otten (1987b), producing the four line segments of the sliding filament theory (dashed line in Fig. A1). The segments were fitted by the formula mentioned above, producing values for a, b, and s (Table A1), with the constraint that their maxima would coincide (Fig. A1). In this way, we could use a single equation for the force-length relationship and could introduce some Gaussian spread in sarcomere length in parallel muscle fibers in a single muscle. The values obtained were: a = 2.46, b = 0.545, and s = 0.216.

Morphometric parameters The morphometric parameters used in our simulations were derived from a data set of one of us (van Eijden). Data were used on sarcomere length, muscle fiber length, and physiological cross-sectional areas as well as three-dimensional coordinates of the attachment sides of three jaw-closing muscles [temporalis (anterior and posterior part), masseter (superficial and deep portion), and medial pterygoid muscle] and the four opening muscles [lateral pterygoid, digastric (anterior and posterior part), mylohyoid and geniohyoid muscle]. The morphometric measurements were taken at a mouth opening of ~3° (equivalent to ~8.5 mm mouth opening including 4 mm overbite; hereafter reference position). The detailed material and methods are described elsewhere (Van Eijden et al. 1995).

We dimensioned the modeled muscles by taking into account the geometric arrangement of the actual jaw muscles. For this, we transformed the data of the human cadaver study of Van Eijden. This was done in the following way.

1)A value for tendon length was obtained by subtracting muscle fiber length from muscle length. (We neglected the angle of pennation of the fiber bundles to the aponeuroses, introducing an error of 3.4% at an estimated angle of pennation of 15°.)

2)We projected the three-dimensional points of attachment of the muscles onto the sagittal plane (Table A2). (By doing so, we introduced an error of 4.4% in element length of the jaw-closing muscles and ~13.8% in the jaw-opening muscles).

3)The resultant maximal bite force (Fbmax) at the level of the canines was derived from
2⋅<LIM><OP>∑</OP></LIM><IT>F</IT><SUB>i⋅</SUB><IT>d</IT><SUB>i</SUB>/<IT>d</IT><SUB>b</SUB> = Fb<SUB>max</SUB>
where Fi is the maximal force of muscle i; di is the length of the moment arm of the muscle vector of muscle i; and db is the length of the moment arm of the bite force vector. This was done separately for the closing muscles and the opening muscles (resulting in Fb1max and Fb2max, respectively). [Maximal forces of the jaw muscles were determined on the basis of their cross-sectional area; we used a value of 0.35 N/mm2 (Nygaard et al. 1983; Weijs and Hillen 1985)]. The moment arms are the perpendicular distances between the two-dimensional force vectors and the axis of rotation of the mandibular condyles.

Because the condyles translate 16 mm anteriorly and 4 mm caudally at a jaw rotation of 30° (Falkenström 1993; Merlini and Palla 1988; Obwegeser et al. 1987), the average position of the axis of rotation was determined to be situated 31.6 mm below and 1.9 mm anterior to the center of the condyles. The length of db was calculated to be 79.2 mm. Fb1max appeared to be 622 N, whereas Fb2max amounted to 288 N) (Table A3).

4)We corrected for the geometric arrangement of the jaw muscles by multiplying the length of the muscles, fibers and tendons by db/di. (This magnifies these lengths by a factor of ~3 for the jaw-closing muscles and a factor of ~2 for the jaw-opening muscles.)

5)The jaw-closing and -opening muscles were grouped by taking the weighted average of their muscle, fiber, and tendon length per group. The weight factors used were the maximal force contribution of each muscle after correcting for the force transmission di/db.

6)From the sarcomere length of each jaw muscle at the jaw reference position, we took a rough average (up to the second decimal) resulting in 2.40 µm for muscle 1 and 2.75 µm for muscle 2. Because the human optimal sarcomere length was found to be 2.7 µm (Fig. A1), we could establish the normalized length of muscles 1 and 2 (0.89 and 1.02, respectively) at the reference position. From this we established that the optimal fiber length of muscle 1 is 75.6 mm and of muscle 2 is 45.1 mm (Table A3).

7)To be able to simulate the covered range of mouth opening of 24.5-33.5 mm, we moved the hyoid bone 5 mm dorsocaudally (without changing the moment arms of the opening muscles) as compared with the position of the hyoid in the material of van Eijden. This choice was based on Pancherz et al. (1986), who showed a similar hyoid movement when the mouth is opened over a distance of 19 mm. After the above geometric transformation, we arrived at a normalized length of muscle 1 of 1.22 and 1.10 and of muscle 2 of 0.71 and 0.91 at a mouth opening of 33.5 and 24.5 mm, respectively (see Fig. 8).

The effective mass of the mandible and the jaw muscles at the level of the incisors was calculated to be 0.230 kg. By effective we mean that we are looking for a mass situated at the incisors that mechanically has the same effect on the dynamics of the motion as the distributed mass of the mandible and jaw muscles. First we calculated the inertia of the mandible and jaw muscles about the axis of rotation in our study from the moment of inertia of the lower jaw and jaw muscles as measured by Koolstra and Van Eijden (1995). When properly chosen, a mass situated at the incisors multiplied by the square distance to the axis of rotation of the lower jaw produces the same inertia as the one calculated.

The values for the internal masses of muscles 1 and 2 in the model were calculated using the total volumes of their muscle fibers and their distances to the center of rotation of the lower jaw (we were interested in the effective masses of the muscles at the bite point). To obtain the volumes, we multiplied the physiological cross-sectional area with the fiber length of each muscle. The masses of the jaw-closing and -opening muscles were calculated to be 0.060 and 0.015 kg, respectively. We divided the masses of the muscles in three equal parts of which one part was linked to the lower jaw, one part included in the model as internal muscle masses (Table A3), and one part was assigned to the unmoving skull or hyoid bone.

Tuned parameters To obtain full definition of the model, a choice still had to be made of another seven parameter values. These parameters describe the passive force-length properties of the muscle fibers and tendinous sheets, the strain of the tendinous sheets at maximal force, the fiber type composition of the muscles and the force production of muscle 2 (Table A4A).

The values of these parameters were calculated by the optimization procedure while meeting two criteria: the model responses should match as close as possible with the experimental results of Nagashima et al. (1997) and the values of the parameters should be within the physiological range as far as this is known. To establish the parameter values, we used the search procedure described by Nelder and Mead (1965). This procedure finds minima of multidimensional functions of which partial derivatives are hard to calculate. Minimizing such functions is impossible with classical steepest descent methods (Press et al. 1986). (Table A4A gives the results of this search.)

Expressions for the passive force-length relationships of both the muscle fibers and tendinous sheets were used as suggested by Otten (1987a).

Identical passive force-length relationships of the muscle fibers were tuned for muscles 1 and 2 after normalizing for both optimallength and maximal isometric active muscle force. After tuning, they gave a passive force of the jaw-closing muscles of ~20 N at a mouth opening of 33.5 mm, which, by coincidence, appeared to be in line with the results of Miles et al. (1986).

Tendinous sheet strain was kept within boundaries of 4-8% while optimizing and was left by the procedure at a value of 5% of their rest length at maximal isometric muscle force. Experimentally, aponeurotic sheet strain is between 6 and 8% of its rest length (Lieber et al. 1991; Trestik and Lieber 1993). This still leaves us with a choice for the stiffness of the tendinous sheets. This value was left free in the fitting process but were kept identical for the sheets of muscles 1 and 2 relative to the maximal forces of the muscles.

In the optimized configuration, the fiber type composition appeared to be 88% type 1 and 12% type 2 for muscle 1 (after applying boundaries of 0-30% type 2), and 17% type 1 and 83% type 2 for muscle 2 (after applying boundaries of 60-90% type 2). The applied boundaries are those described by Eriksson et al. (1981, 1982) and by Eriksson and Thornell (1983).

As a result of the optimization procedure, the force exerted by the cocontracting jaw-opening muscles appeared to be 5.1 N.

Damping properties were added to the tendinous sheets because hysteresis is a well known property of tendinous material and because they were convenient and necessary in avoiding intramuscular oscillations in the simulations. We chose a low damping coefficient of the tendinous sheets of 5 N s/m (Table A4B).

Device parameters The contribution of each of the four forces (Foffset, Fmagnet, Fvibration and Fend) to Fresistance and the effective mass of the lower bar at the bite point was established by the above-mentioned search method in which the behavior of the device model was matched with that of the unloading apparatus. The effective mass of the lower bar at the bite point was 1.40 kg.

    FOOTNOTES

  Address for reprint requests: G.E.C. Slager, Dept. of Medical Physiology, University of Groningen, Bloemsingel 10, 9712 KZ Groningen, the Netherlands.

  Received 15 April 1997; accepted in final form 31 July 1997.

    REFERENCES
Abstract
Introduction
Methods
Results
Discussion
References

0022-3077/97 $5.00 Copyright ©1997 The American Physiological Society