Determination of the in situ mechanical behavior of ankle ligaments

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

Available online at www.sciencedirect.com

www.elsevier.com/locate/jmbbm

Research Paper

Determination of the in situ mechanical behavior of ankle ligaments Bingbing Nien, Matthew B. Panzer, Adwait Mane, Alexander R. Mait, John-Paul Donlon, Jason L. Forman, Richard W. Kent Center for Applied Biomechanics, University of Virginia, 4040 Lewis and Clark Drive, Charlottesville, VA 22911, USA

ar t ic l e in f o

abs tra ct

Article history:

The mechanical behavior of ankle ligaments at the structural level can be characterized by

Received 15 March 2016

force–displacement curves in the physiologic phase up to the initiation of failure. However,

Received in revised form

these properties are difficult to characterize in vitro due to the experimental difficulties in

26 August 2016

replicating the complex geometry and non-uniformity of the loading state in situ. This

Accepted 7 September 2016

study used a finite element parametric modeling approach to determine the in situ

Available online 13 September 2016

mechanical behavior of ankle ligaments at neutral foot position for a mid-sized adult foot

Keywords:

from experimental derived bony kinematics. Nine major ankle ligaments were represented

Ankle ligaments

as a group of fibers, with the force–elongation behavior of each fiber element characterized

Mechanical behavior

by a zero-force region and a region of constant stiffness. The zero-force region, represent-

In situ Parametric finite element modeling Collagen fibers Biomechanics

ing the initial tension or slackness of the whole ligament and the progressive fiber uncrimping, was identified against a series of quasi-static experiments of single foot motion using simultaneous optimization. A range of 0.33–3.84 mm of the zero-force region was obtained, accounting for a relative length of 6.773.9%. The posterior ligaments generally exhibit high-stiffness in the loading region. Following this, the ankle model implemented with in situ ligament behavior was evaluated in response to multiple loading conditions and proved capable of predicting the bony kinematics accurately in comparison to the cadaveric response. Overall, the parametric ligament modeling demonstrated the feasibility of linking the gross structural behavior and the underlying bone and ligament mechanics that generate them. Determination of the in situ mechanical properties of ankle ligaments provides a better understanding of the nonlinear nature of the ankle joint. Applications of this knowledge include functional ankle joint mechanics and injury biomechanics. & 2016 Elsevier Ltd. All rights reserved.

n

Corresponding author. Fax: þ1 434 297 8083. E-mail address: [email protected] (B. Nie).

http://dx.doi.org/10.1016/j.jmbbm.2016.09.010 1751-6161/& 2016 Elsevier Ltd. All rights reserved.

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

1.

Introduction

Syndesmosis and lateral ligament sprains account for about 85–90% of foot and ankle injuries (Rubin and Sallis, 1996), which are considered the most vulnerable regions of the lower extremity in all levels of sports (Barker et al., 1997; Bloemers and Bakker, 2006; Hinterman, 1999). Stability of the ankle joint is dictated by the ligaments and tendons as a restraint, emphasizing the significance of a precise identification of the structural properties of ankle ligaments. Such properties are referred to as the “in situ” behavior since they represent the ligament responses in the normal anatomic position in the ankle. The literature still lacks a precise characterization on the mechanical behavior of ankle ligaments in situ under gross foot and ankle kinematics. At the structural level, a ligament can either be under tension (i.e. pretension) or slack (i.e. resulting in a zero-tensile force) for a given joint position (Fleming and Beynnon, 2004). For example, Ozeki et al. (2002) reported that lateral ligaments of the ankle were slack around the neutral foot position (the right-angle position in the tibio-pedal angle), exhibiting an initial elongation of 3% to  2% compared to the zero-force length. Nigg et al. (1990) estimated that the initial slackness of some deltoid ligaments would account for 40% of the ultimate strain. At material level, ankle ligaments are soft connective tissues consisting of densely packed collagen fibers that preferentially aligned along the longitudinal axis to transfer load among bones (Leardini et al., 2000). Collagen fibers are crimped to varying degrees while unloaded and a fiber's crimping is progressively straighten as load is applied, accounting for a low-force region, known as the toe region. The loading region, which follows the toe region and can be idealized as linear (i.e., as having constant stiffness) up to the initiation of fiber failure, corresponds to the stretching of uncrimped collagen fibers to build ligament force rapidly with progressive lengthening (Lucas, 2008; Woo et al., 1999). The potential interactions among fibers are reported to be mechanically insignificant at low strain rates (Piérard and Lapière, 1987). Therefore, for a given joint position, the lump mechanical behavior of individual fiber element can be characterized by nonlinear force–displacement curves with distinct zero-force and loading regions up to the initiation of failure (Decraemer et al., 1980; Lucas, 2008; Woo et al., 1999). The zero-force region represents this initial tension or slackness at the structural level and the progressive uncrimping as an inherent material property. For the purpose of biomechanical analysis, individual ligament was composed of a collection of fibers with similar properties. In the physiologic phase, the ligament maintains structural integrity and function normally without damage, i.e., controls joint kinematics at the gross level (Mattucci and Cronin, 2015; Yoganandan et al., 2001). Numerous research efforts have been made to understand the biomechanical response of ankle ligaments. Experimental studies are commonly performed using isolated boneligament-bone specimens to characterize ligament tensile behavior (Funk et al., 2000; Hall, 1998; Johnson and Markolf, 1983; Rasmussen, 1985; Xenos et al., 1995). Recent studies using robotic technology have offered a practical way of

503

determining in situ forces of individual ligament subjected to pre-defined joint kinematics, particularly for knee ligaments (Fox et al., 1998; Gabriel et al., 2004; Takai et al., 1993; Zhang et al., 2016). These studies were generally based on the principle of superposition, i.e., the difference between the force data recorded in the intact joint and the joint with some dissected ligament was believed to be the in situ force of the removed ligament. The potential influence of removing supporting structures on joint mechanics has not been taken into account. Overall, the difficulty of measuring the inhomogeneous force–elongation distributions and load paths of multiple ligaments has precluded quantitative description of in situ ligament behavior that provides joint stability within an intact joint (Funk et al., 2000; Mommersteeg et al., 1997). In addition to experimental studies, finite element (FE) models provide the capability to investigate the stress and strain built in soft tissues and are expected to replicate the biomechanical responses of the ankle (Cheung et al., 2006; Reggiani et al., 2006). When existing modeling implementations have been used as a flexible way to study ankle kinematics (Tannous et al., 1996; Wei et al., 2011), the ligamentous structure and complex loading path have not been sufficiently considered yet. Comprehensive description on the structural response of ankle ligaments is needed to characterize the in situ loads and the associated elongation along the longitudinal axis. Therefore, the objective of this study was to determine the in situ mechanical behavior built in the microstructure of ankle ligaments at physiological levels of loading. A previous developed parametric framework for FE modeling of ankle ligaments was utilized as a link between gross motions of the ankle and subtalar joints and the underlying ligament mechanics that generate them. The zero-force region and the stiffness within the loading region of ankle ligaments were identified via optimization techniques in minimizing the difference between the ankle model and experimental derived kinematics. The ankle model with implemented in situ ligament behavior was then evaluated based on the kinetic and kinematic responses under multiple loading conditions. It is anticipated that comprehensive replications of gross foot motions can, by implication then, identify the in situ deformations of the ligaments at the microstructural level.

2.

Materials and methods

We used a simultaneous optimization approach to investigate the contribution of each ligament and their combinations systematically along with the computational foot and ankle model and the experimental derived bony kinematics (Fig. 1). An existing parametric model was utilized as the basis of optimization (Nie et al., 2016). Four independent suboptimizations on different single rotations of the foot were performed simultaneously. For each suboptimization, the simulation was setup replicating the experimental boundary constraints. The obtained optimum, representing the in situ ligament behavior, was evaluated considering the predicted ligament responses under multiple foot rotations. A detailed description of the computational framework, the

504

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

experimental design, the simulation setup and the model evaluation is presented in the following subsections.

The nine major ligaments were numbered from 1 to 9 (i¼ 1, …, 9) and each was represented as a group of fibers. Fibers were modeled by one-dimensional (1-D) discrete ele-

2.1.

Computational framework

The foot and ankle model used in this study represents a typical adult foot, with the bone geometry precisely reconstructed from a mid-sized male volunteer (175.3 cm height and 77.1 kg weight). The kinematic response of the model was previously validated under different loading conditions, including quasi-static axial rotations (Shin et al., 2012). The bones were treated as rigid bodies. Modeling of deformable articular cartilage was attached to bones to reduce the gaps at bony interface and to provide physically realistic joint kinematics. A highly compressible elastic material model was defined, with an elastic modulus of 12 MPa, a nominal density value of 1000 kg /m3 and a nonlinear stress–strain curve adopted from (Mow and Mansour, 1977). A kinematic constraint-based contact algorithm was implemented for cartilage-to-cartilage contact and a penalty-based surfaceto-surface contact algorithm was assigned to allow for compression loads to be transferred among bones. Given the fact that joint kinematics is mainly guided by ligaments and the articulating geometry (Leardini et al., 2000; Weiss and Gardiner, 2001), flesh modeling was excluded in the present ankle model to improve computational efficiency. Nine major ligaments were replicated as a bundle of collagen fibers: the anterior talofibular (ATaF), posterior talofibular (PTaF) and calcaneofibular (CF) ligaments on the lateral side; the anterior tibiotalar (ATaT), posterior tibiotalar (PTaT), tibiocalcaneal (CT), and tibionavicular (TiN) ligaments, which composed the deltoid on the medial side; and the anterior tibiofibular (ATiF) and posterior tibiofibular (PTiF) ligaments that comprise the distal syndesmosis.

ments, with two ends rigidly attached to ankle bones, and evenly distributed along the pre-determined insertion widths based on dissection studies and human anatomy (Fig. 2(a)). The overall mechanical behavior of each fiber element was characterized by a bilinear force-displacement curve with two zones: a zero-force region, which was defined in unit of length, and a region of constant stiffness (Fig. 2(b)). Within each ligament, the bundle of fiber elements had the same initial zero-force region in absolute terms, which was defined as a variable, ci (mm) (i¼1, …, 9). Ligaments are made up of varying amounts of the same microstructural components. Therefore, the number of the fibers included in each ligament was recorded as ni , with the same stiffness of 0.02 kN/mm adopted across all the ligaments (Lucas, 2008). As a result, the total stiffness of each ligament under uniaxial tension, Ki, can be estimated as Ki ¼ ni *0.02. In total, eighteen variables, ci and ni (i¼ 1, …, 9), were included in the parametric model characterizing the ligamentous behavior at the microstructural level. The variables were recorded as an s-vector X¼ (x1, …, xs) (s¼ 18). For a given value of X, an ankle model could be generated using scripts written in MATLAB (R2015b, The Math Works Inc., Natick, MA, USA) with the specified zero-force region and fiber number for each ligament. The value of each variable, xp (p ¼1, …, s), was identified by matching the structural response predicted by the model with the experimental data through optimization (Fig. 1). The overall optimum results of X were obtained as the in situ ligament behavior. The gross structural response of ankle ligaments was presented by the total force–displacement responses. The total tension force indicated the summation of force sustained by all the fibers within each ligament; the

Fig. 1 – Flowchart of the optimization for determining in situ ligament behavior.

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

Fig. 2 – (a) Parametric modeling of the ligament in the foot and ankle model (shown: the ATiF, ATaF and CF ligaments). (b) Overall mechanical behavior of one fiber element within each ankle ligament. elongation indicated the average elongation of the fiber elements.

2.2.

Experimental data overview

Bony kinematics under gross foot rotation was obtained using the experimental design previously reported by Mait et al. (2015). The lower extremities from three human cadavers (mean age, 42 years; range, 31–47 years) were disarticulated at the knee taking care to retain the proximal ligaments between the tibia and fibula. The proximal portion of the tibia was potted in a cylindrical cup using wood screws and Bondo body filler (Part #261, 3M Company, St. Paul, MN, USA). The motion of the fibula was dictated only by its connective tissue (i.e., the fibula was not potted). The calcaneus was potted and installed in a jig that allowed it to translate in the transverse plane (i.e., along the X- and Y-axis, Fig. 3a). The foot was initially positioned in a neutral orientation, i.e., the distal end of the first phalanx, the centroid of the calcaneus, and the long axis of the tibia approximately formed a 90-degree angle in the sagittal plane. No axial loading was applied to the initially positioned foot. To produce a single rotation of the foot, rotation loading was applied to the calcaneus quasistatically and the corresponding applied moment was measured by a load cell mounted on top of the potted tibia. Then

505

Fig. 3 – (a) The experimental test rig and the measurement of bony kinematics. (b) Illustration of the Local Coordinate System (LCS) defined for each bone (shown: left foot).

Table 1 – Experimental loading conditions and the input rotation angles of the calcaneus (subject 743L). Loading conditions

Rotation angle (deg)

Single rotation (for model optimization) External rotation Internal rotation Eversion Dorsiflexion

20.7 25.38 11 16.81

Multiple loadings (for model evaluation) Eversionþexternal rotation Dorsiflexionþexternal rotation

11.00þ20.50 16.81þ22.55

the loads were removed and the foot returned to its neutral orientation. Four different single rotations were conducted, including external rotation, internal rotation, eversion and dorsiflexion. Multiple loadings were performed by applying different single rotations in sequence (e.g. an eversion followed with an external rotation) and were used for model evaluation purpose (Table 1). To provide a quantitative description of the relative motions among ankle bones, we used Local Coordinate Systems (LCS) for the ankle bones including the tibia, fibula, talus, calcaneus and navicular with the procedure proposed

506

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

by Camacho et al. (Fig. 3a). Such approach proved capable in capturing the subtle motion among the relative small bones, e.g. the navicular (Camacho et al., 2002). The LCSs were defined using aluminum cubes mounted via screws to each of the bones. During the test, each cube was scanned at the initial position and final position associated with the rotation of the foot using a 7-axis, three-dimensional laser scanner (ROMER Absolute Arm, Hexagon Metrology, Inc., Surrey, Great Britain). Using the scanned data, the tibia coordinate system served as the reference frame; the resultant changes in yaw– pitch–roll angle, i.e., the Euler angle description (z–y–x), and the x, y, and z displacements (recorded as Δdx, Δdy and Δdz) of the fibula, talus, navicular and calcaneus, relative to the tibia, was calculated from the direction cosine matrix (Fig. 3b). The subject 743L was anatomically closest to the model among the three specimens, so it was chosen as the experimental basis for the model ligament optimization. The input rotation angles to the calcaneus of subject 743L under the single rotations and multiple loadings were provided in Table 1. Under each rotation input from the calcaneus, the difference on bony kinematics predicted by the model and measured from the experiment was quantified by KINdiff, which was defined as the mean value of two normalized root mean squared error of the nine bony angles and displacements, as below ffi 0sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 9  P 1 Angle Angle B 9 model;r experiment;r B r¼1 KINdiff ¼ B B maxjAngleexperiment;r j @ sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1 9  2 P 1 Δdmodel;r  Δdexperiment;r C 9 C r¼1 C=2 þ C maxjΔdexperiment;r j A

ð1Þ

where r is the number of the bony angle or bony displacements, i.e. r ¼1, …, 9. The subscript “model” and “experiment” represented the computationally and experimentally   obtained bony angles  or displacements; maxAngleexperiment;r    and maxΔdexperiment;r  denoted the maximum value of the experimental data. Given the difference in morphometry of the model and the three specimens, the minimum and maximum values of the bony angles and displacements of all three subjects were considered to be an indication of the experimental range and thus provide context for the model's response.

2.3.

Simulation setup

The experimental setup was replicated using the FE ankle model in explicit LS-DYNA software (mpp971sR7, LSTC, Livermore, CA, USA). The proximal end of the tibia was fixed in space, and the calcaneus was positioned to match each experimental condition. Continuous rotation was applied to the calcaneus to enforce a prescribed rotation at 101/s to the angle measured in the experiments (Fig. 3b, Table 1). The fibula, talus and navicular were free in all six DOF so that their motion was dictated by the relevant connective tissues, including ligaments. The resultant yaw–pitch–roll angles of the fibula, talus and navicular were determined from

coordinate transformation and consequently, a total of nine angles (three angles x three bones) and nine displacements (x, y, and z displacements x three bones), defined analogously to the experiments, were available as the model response. Within each ligament, the average length of the bundle of fiber elements was used to represent the initial length, l0,i. The variable range (upper and lower limits) of the zero-force region (ci) in the optimization was set as 0–20% of the initial ligament length, l0,i. The initial stiffness of each ligament, K0,i, was estimated via linear approximation of the rapid increase in force generation along the force–displacement behavior as reported in the tensile tests of bone-ligament-bone specimens (Funk et al., 2000). This resulted in the initial value of the fiber number being determined as n0,i ¼ K0,i/0.02. The variable range of ni was set as 0–2 times of the initial value of n0,i. Four suboptimizations, corresponding to the four different single foot rotations in the experiments, were recorded as sopt (m) (m¼1, …, 4): (1) external rotation; (2) internal rotation; (3) eversion and (4) dorsiflexion (Fig. 1). The suboptimizations were performed simultaneously to reproduce the bony kinematics. The Kriging-based genetic algorithm, which yielded high-potential to create accurate global approximations to facilitate multi-dimensional design optimization, was used (Simpson et al., 2001). The quantitative kinematics difference between the model and the experiment, KINdiff, was used as the objective function to minimize in the optimization. The optimization was performed using LS-OPT (V5.1.0, LSTC, Livermore, CA, USA). After each iteration, convergence of the variables from each suboptimization, xp, sopt(m) (p¼ 1, …, s, m¼ 1, …, 4), was assessed and those converged variables were set as constant in the next iteration. The iterations were stopped when the values of the objective function converged, i.e., the norm of the gradient of the objective function was lower than 0.0001.

2.4. Evaluation of the model implemented with in situ ligament behavior To evaluate the ankle model implemented with the optimal ligament behavior, simulations under the two sequential multiple loading conditions were conducted to compare with the experiment (Table 1). One included an 11.001 eversion followed by a 20.501 external rotation of the foot. The other included 16.811 of dorsiflexion and 22.551 of external rotation applied in a similar manner. The bony angle and displacements of the fibula, talus and navicular were measured. The applied rotation moment to produce the prescribed calcaneus motion was compared with the experimental range. Comparison of the bony kinematics of the experimental and simulation results was quantified using the pre-defined indicator, KINdiff, to assess the prediction capability of the ligament modeling in the loading conditions other than those used in the optimization.

3.

Results

3.1.

In situ mechanical behavior

Convergence of the objective function, KINdiff, in the four suboptimizations was achieved after about 140 iterations

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

(Fig. 4). Zero-force regions ranging from 1.4871.22 mm were obtained from the optimal results, which corresponded to

507

gross structural response of the ligaments, i.e., the force– displacement responses, was presented in Fig. 6.

relative lengths of 6.773.9% of the initial ligament length (Table 2). On the medial side, the PTaT contributed a sig-

3.2.

Model evaluation

nificant portion to the deltoid stiffness and exhibited a zeroforce region of 1.73 mm, while the TiN ligament tended to provide a low-stiffness of 0.039 kN/mm with a zero-force region being up to 3.84 mm. The two ligaments at the tibiofibular syndesmosis, ATiF and PTiF, exhibited a zeroforce region of 1.25 mm and 0.38 mm respectively, followed with a high-stiffness. Bony kinematics of the ankle model with the optimum variable combination was generally in accordance with the experimental data, with a kinematics difference (KINdiff) of 0.154, 0.155, 0.182, 0.150 under external rotation, internal rotation, eversion and dorsiflexion of the foot, respectively (Fig. 5). Some deviations were noticed on navicular roll under internal rotation (model, 6.431 experimental range,  19.00 to 21.491 Subject 743L,  19.731) and navicular yaw under eversion (model,

1.801 experimental range,

0.39 to

 12.681 Subject 743L,  6.491), with large variability also being

In response to the multiple loading conditions, the computational model exhibited a 3.0% difference on the applied moment (28.7 Nm vs. 29.6 Nm for Subject 743L) to generate the prescribed calcaneus motion when compared to the experiment in the “eversionþexternal rotation” loading condition. The peak moment predicted by the model (31.9 Nm) was 1.3 Nm beyond the experimental range (16.2–30.6 Nm) in the “dorsiflexionþexternal rotation” loading condition (Fig. 7). A kinematics difference (KINdiff) of 0.116 and 0.153 was obtained when compared to subject 743L in the “eversionþexternal rotation” and “dorsiflexionþexternal rotation” loading condition (Fig. 8). Most of the bony angles and CG displacement predicted by the model were within the experimental range. This indicated a good biofidelity of the ankle model implemented with in situ ligament behavior under multiple loading conditions other than the single rotational loading used in the optimization.

observed among the experimental results. Following this, the

Fig. 4 – Histories of the objective function in the four suboptimizations on different foot rotations.

4.

Discussion

4.1.

Determination of the in situ ligament behavior

This study determined the in situ nonlinear mechanical behavior of ankle ligaments with considerations of the zero-force region and loading region as an inverse-problem of parameter identification. The ligament parameters were obtained via optimization in an effort to minimize the difference between one representative ankle model and experimental derived kinematics. The optimization converged to the solution within the defined searching range. The novel computational approach used in this study demonstrates the feasibility to explore the in situ state within the ligament and exhibits some unique

Table 2 – In situ behavior of ankle ligaments. Ligament list

No. (i)

Ligament geometry

Zero-force region

Loading region

Insertion width (mm)

l0,i (mm)

ci (mm)

ci (%)

ni

Ki (kN/mm)

Deltoid (Medial)

ATaT PTaT CT TiN

1 2 3 4

6.4 (Talus), 7.0 (Tibia) 6.4 (Talus), 7.0 (Tibia) 6.0 (Calcaneus), 6.0 (Tibia) 12.6 (Tibia), 9.3 (Navicular)

20.8 14.3 24.3 33.2

0.97 1.73 2.34 3.84

4.7 12.1 9.6 11.6

6 41 15 2

0.125 0.821 0.302 0.039

Lateral

ATaF PTaF CF

5 6 7

3.5 (Talus), 4.7 (Fibula) 4.8 (Talus), 4.3 (Fibula) 10.0 (Calcaneus), 10.0 (Fibula)

16.5 25.0 21.4

0.33 1.65 0.68

2.0 6.6 3.2

4 39 23

0.089 0.775 0.452

High ankle

ATiF PTiF

8 9

15.0 (Tibia), 15.0 (Fibula) 15.0 (Tibia), 10.0 (Fibula)

16.5 14.8

1.25 0.38

7.6 2.6

43 47

0.870 0.934

Note: l0,i (mm) indicated the initial length of the ligaments; ci (mm) indicated the optimum value of the zero-force region; ci (%) indicated the relative length of the optimum value, which was calculated as its ratio relative to the initial length of the ligament. Ki indicated the estimated linear stiffness of the ligament under uniaxial tension.

508

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

Fig. 5 – Comparison of the bony angles of the experiments and the model on single rotation of the foot (model optimization).

Fig. 6 – Gross structural response of ankle ligaments.

features. Without the necessity of external devices being attached to the soft tissues in the experimental setup, it eliminated the potential risk of altering the load–elongation behavior of the ligaments. An ankle model validated at the ligament level, such as the one built up in this study, can

characterize the effective non-uniform force–displacement states that are not accessible to experimental observations. Modeling the individual ligament as a bundle of 1-D fiber elements is a reasonable balance of the gross mechanical behavior of ligaments, the microstructure responsible for it, and the bony kinematics that result. Each fiber's mechanical behavior was characterized using an idealized bilinear curve (Fig. 2(b)). Such an approach also provides a framework of other constitutive models of varying complexity for collagen fibers, such as viscoelastic models. The gross structural response of ankle ligaments derived from the fiber elements (Fig. 6(a)) was in accordance with available data in existing experimental studies. For the CF and PTaF ligaments on the lateral side, the total tension force reached approximately 118 N and 177 N with 1 mm average elongation in the loading region. Comparably, Siegler et al. (1988) estimated a linear stiffness of 126.6 N/mm and 164.3 N/mm after toe region for these two ligaments. For the ATiF ligament, our model had a zero-force region of 7.6%. In the quasilinear viscoelasticity model of the ATiF ligament proposed by Funk et al. (2000), a noticeable force increase was after an engineering strain of about 6%. Some deviations were noticed on the navicular rotation angle when optimizing the ligament behavior under single foot rotation. Among the investigated ankle ligaments, only the TiN ligament was attached to the navicular.

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

509

The navicular is a small bone located on the medial side of the foot, and articulates proximally with the talus. The TiN ligament was observed as a thin, homogeneous band of fibers attached from the tibia to the navicular and inserted posteriorly (Muhle et al., 1999). The low-stiffness and high zeroforce region of the TiN ligament (Table 2) led to difficulties in replicating the navicular kinematics accurately. Nevertheless, the ankle model implemented with in situ ligament behavior proved capable of predicting the kinematics response under multiple loading conditions (Fig. 8). The trends in ligament behavior observed with respect to the bony angle and translation are believed to be overall reasonable.

4.2.

Fig. 7 – Comparison of the applied rotation moment of the experiments and the model under multiple loading conditions.

In situ mechanical behavior of ankle ligaments

Knowledge of the in situ mechanical behavior of ankle ligaments is vital to explain their role in ankle joint function. The sequential uncrimping of the collagen fascicles during ligament tension before the development of tensile force has estimated to occur at a strain range of 2–7% (Viidik, 1990), depending on the collagen content and morphology of individual ligaments. Besides, the in situ state and effective

Fig. 8 – Comparison of bony kinematics of the experiments and the response predicted by the model on multiple loading conditions (model evaluation).

510

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

loading paths at the structural level of a foot complex are difficult or impossible to measure experimentally (Weiss and Gardiner, 2001). This uncertainty complicates the study of ankle mechanics and was the primary motivation of the current study. The sequential uncrimping of fibers and the initial in situ states of ligaments were lumped in to the zeroforce region in the present study, and proved to play an important role in dictating bony kinematics. The anterior superficial deltoid, i.e., the TiN ligament, exhibited a zero-force region range of 11.6% and relative low-stiffness compared to the posterior deltoid ligaments (Table 2), which implies a late involvement to resist excessive rotation at the subtalar joint. For example, the normal range of plantar flexion was reported to be 10–551 (mean 39.61), which is higher than dorsiflexion (mean 15.31) (Roaas and Andersson, 1982). The TiN ligament helps to restrict plantar flexion of the foot, where the talus moves progressively forwards out of the mortise by the malleoli and the joint. Such behavior of the TiN ligament contributes to the highdegree of freedom in addition to the shape of the bony surface. The lateral ankle ligaments are at a high-risk of sport injuries (Anderson et al., 2010), with the ATaF ligament as the most commonly injured ligament in lateral ankle sprains, followed by the CF ligament (Ray et al., 1997). Existing anatomical dissections suggested that ATaF and CF ligaments contributed a significant portion to the lateral ankle stability (Heilman et al., 1990; Michael et al., 2008; Weindel et al., 2010). Such function is influenced by a number of important parameters including mechanical properties, initial strain and anatomical geometry. The ATaF ligament exhibited a short zero-force region of 0.33 mm (2.0% of the initial length), indicating an early increase of tensile force within the ATaF ligament during ankle motions that impose load on it, e.g., the restraint that ligament applies to the talus under flexion or internal rotation of the foot at the ankle joint. The “high” ankle ligaments (ATiF and PTiF) where characterized as having relatively short zero-force regions of 7.6% and 2.6% and high stiffness (Table 2). This supports the functional role of these two ligaments as the primary constraint to excessive fibular motion relative to the tibia near the syndesmosis (Lin et al., 2006; Williams et al., 2007). This may be an explanation for the clinical observation that syndesmosis injuries often occur in conjunction with deltoid ligament injury and/or fractures in the malleoli (Fallat et al., 1998; Hopkinson et al., 1990). Xenos et al. (1995) reported that sectioning of the ATiF ligament resulted in diastasis of 2.3 mm at the distal tibiofibular joint, which is consistent with the 1.25 mm zero-force region followed by a rapid force increase found in this study.

4.3.

forces have not been considered. In addition to the idealized bi-linear fiber behavior used in the present study, other phenomenological or continuum material models, such as the quasilinear viscoelasticity model commonly applied in biomechanics (Forestiero et al., 2014; Fung, 1972; Funk et al., 2000; Hall, 1998; Subit et al., 2008), can be used in subsequent studies to explore additional aspects of ligament behavior during various loading conditions. The bone geometry was reconstructed from a mid-sized male volunteer and ligament insertions were anatomically reasonable but nominal. The model is assumed as a representation of a typical adult foot; the experimental data examined in the optimization efforts were limited to the selected specimens. Biological variability was not accounted for and the degree to which it contributes to the noticed deviations between the optimal model and the experimental results (Fig. 5) remains unknown. Therefore, a preliminary parametric study was performed on the ATiF ligament to analyze the sensitivity of ligament behavior to the selected geometry (e.g. insertion area) and the initial ligament strain (e.g. the zero-force region) obtained from optimization. The insertion width and the optimal zero-force region of the ATiF ligament (referred by as c(ATiF)) were varied by 720%, respectively. Minimal effects were found when changing the insertion width by 20%. The resultant change in ligament behavior due to varying c(ATiF) were compared in Fig. 9. When reaching a total tension force of 400 N, varying c(ATiF) by 10% or 20% led to a difference of approximate 0.1–0.2 mm on average ligament elongation. The trends of the tension force and elongation of the ligament with the given bony kinematics were similar. Overall, the changes were not significant, indicating that the reported optimal ligament parameters are not likely to be dramatically affected by small variations on the ligament geometry or initial strain. Given the know effects of age, sex and anthropometry on the absolute value of the in situ behavior in the ankle ligaments, further research efforts are necessary to quantify the potential influence on the ligament response due to the intersubject variability in morphometry to a broad range. Finally, the ligament parameters at the fiber level were identified via optimization with an effort to minimize the kinematic difference between the representative ankle model and the experiments. It needs to be noted that optimization with high non-linearity is usually not a convex problem, i.e., there can be more than one optimal solution subjected to the optimization technique and design space. Therefore, the reported optimal ligament parameters, which allowed the

Limitations and future research

The feasibility of the method in determining the in situ ligament behavior was demonstrated in the selected loading scenario. The obtained results represent the basic ligament behavior under limited loading rates in the experimental design and therefore must be interpreted with exercise caution. The consequences of higher loading rates, greater rotation magnitudes, and the weight bearing and muscle

Fig. 9 – Gross structural behavior of the ATiF ligament when varying c(ATiF) by 10% and 20%.

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

model to formulate the measured joint kinematics with trustable range, should be considered as a representative but may not be a unique solution. Despite the expected limitations, the data obtained on in situ behavior of the ankle ligaments was useful in better characterizing the function of individual ankle ligaments. This information can help in accurate computational modeling of the ankle to study effective strain in ankle ligaments which are not accessible to experiments. Among the possibilities for further work involving this approach, one of the most significant extensions of the presented study will involve identifying the injury patterns of ankle ligaments at the fiber level.

5.

Conclusions

This paper determined the in situ mechanical behavior built in the microstructure of ankle ligaments under gross foot motion via optimization based on experimental derived kinematics. The parametric modeling of the ligament property demonstrated the feasibility to link the gross structural behavior and the underlying bone and ligament mechanics. The ankle model implemented with in situ ligament behavior proved capable of accurately predicting the bony kinematics under multiple physiological loading conditions. Determination of the in situ mechanical properties of ankle ligaments led to an enhanced understanding of the nonlinear nature of the ankle joint and can ultimately help to elucidate ankle mechanics, identify injury mechanisms, and improve the design of injury countermeasures and orthopedic devices.

Acknowledgments One of the co-authors (Kent) has an ownership interest in one of the sponsors of the study (Biocore). The authors would like to thank the members of the Foot & Ankle Subcommittee of the National Football League for funding, supporting and providing valuable input to this study.

r e f e r e n c e s

Anderson, R.B., Hunt, K.J., McCormick, J.J., 2010. Management of common sports-related injuries about the foot and ankle. J. Am. Acad. Orthop. Surg. 18, 546–556. Barker, H.B., Beynnon, B.D., Renstrom, P.A., 1997. Ankle injury risk factors in sports. Sport. Med. 23, 69–74. Bloemers, F.W., Bakker, F.C., 2006. Acute ankle syndesmosis injury in athletes. Eur. J. Trauma 32, 350–356. Camacho, D.L.A., Ledoux, W.R., Rohr, E.S., Sangeorzan, B.J., Ching, R.P., 2002. A three-dimensional, anatomically detailed foot model: a foundation for a finite element simulation and means of quantifying foot-bone position. J. Rehabil. Res. Dev. 39, 401–410. Cheung, J., An, K., Zhang, M., 2006. Consequences of partial and total plantar fascia release: a finite element study. Foot Ankle Int. 27, 125–132. Decraemer, W.F., Maes, M.A., Vanhuyse, V.J., Vanpeperstraete, P., 1980. A non-linear viscoelastic constitutive equation for soft

511

biological tissues based upon a structural model. J. Biomech. 13, 559–564. Fallat, L., Grimm, D.J., Saracco, J.A., 1998. Sprained ankle syndrome: prevalence and analysis of 639 acute injuries. J. Foot Ankle Surg. 37, 280–285. Fleming, B.C., Beynnon, B.D., 2004. In vivo measurement of ligament/tendon strains and forces: a review. Ann. Biomed. Eng. 32, 318–328. Forestiero, A., Carniel, E.L., Natali, A.N., 2014. Biomechanical behaviour of ankle ligaments: constitutive formulation and numerical modelling. Comput. Methods Biomech. Biomed. Eng. 17, 395–404. Fox, R.J., Harner, C.D., Sakane, M., Carlin, G.J., Woo, S.L.-Y., 1998. Determination of the in situ forces in the human posterior cruciate ligament using robotic technology a cadaveric study. Am. J. Sports Med. 26, 395–401. Fung, Y.C., 1972. Stress strain history relations of soft tissues in simple elongation. In: Fung, Y.C., Perrone, N., Anliker, M. (Eds.), Biomechanics, Its Foundations and Objectives. Prentice Hall, Englewood Cliffs, NJ. Funk, J.R., Hall, G.W., R, C.J., Pilkey, W.D., 2000. Linear and quasilinear viscoelastic characterization of ankle ligaments. J. Biomech. Eng. 122, 15–22. Gabriel, M.T., Wong, E.K., Woo, S.L.Y., Yagi, M., Debski, R.E., 2004. Distribution of in situ forces in the anterior cruciate ligament in response to rotatory loads. J. Orthop. Res. Off. Publ. Orthop. Res. Soc. 22, 85–89. Hall, G.W., 1998. Biomechanical Characterization and Multibody Modeling of the Human Lower Extremity, School of Engineering and Applied Science. University of Virginia, United States. Heilman, A.E., Braly, W.G., Bishop, J.O., Noble, P.C., Tullos, H.S., 1990. An anatomic study of subtalar instability. Foot Ankle Int. 10, 224–228. Hinterman, B., 1999. Biomechanics of the unstable ankle joint and clinical implications. Med. Sci. Sports Exerc. 31, S459–S469. Hopkinson, W.J., Pierre, P.S., Ryan, J.B., Wheeler, J.H., 1990. Syndesmosis sprains of the ankle. Foot Ankle Int. 10, 325–330. Johnson, E.E., Markolf, K.L., 1983. The contribution of the anterior talofibular ligament to ankle laxity. J. Bone Jt. Surg. Am. 65, 81–88. Leardini, A., O’Connor, J.J., Catani, F., Giannini, S., 2000. The role of the passive structures in the mobility and stability of the human ankle joint: a literature review. Foot Ankle Int. 21, 602–615. Lin, C.F., Gross, M.L., Weinhold, P., 2006. Ankle syndesmosis injuries: anatomy, biomechanics, mechanism of injury, and clinical guidelines for diagnosis and intervention. J. Orthop. Sports Phys. Ther. 36, 372–384. Lucas, S.R., 2008. A Microstructural Viscoelastic Ligament Failure Model, School of Engineering and Applied Science and the School of Medicine. University of Virginia, United States. Mait, A., Donlon, J.P., Mane, A., Forman, J., Kent, R., 2015{C}{C}. Kinetics and kinematics of the ankle during foot external rotation. In: Proceedings of the 39th Annual Meeting of the American Society of Biomechanics. Columbus, Ohio, United States. Mattucci, S.F., Cronin, D.S., 2015. A method to characterize average cervical spine ligament response based on raw data sets for implementation into injury biomechanics models. J. Mech. Behav. Biomed. Mater. 41, 251–260. Michael, J.M., Golshani, A., Gargac, S., Goswami, T., 2008. Biomechanics of the ankle joint and clinical outcomes of total ankle replacement. J. Mech. Behav. Biomed. Mater. 1, 276–294. Mommersteeg, T.J.A., Huiskes, R., Blankevoort, L., Kooloos, J.G.M., Kauer, J.M.G., 1997. An inverse dynamics modeling approach to determine the restraining function of human knee ligament bundles. J. Biomech. 30, 139–146.

512

journal of the mechanical behavior of biomedical materials 65 (2017) 502 –512

Mow, V.C., Mansour, J.M., 1977. The nonlinear interaction between cartilage deformation and interstitial fluid flow. J. Biomech. 10, 31–39. Muhle, C., Frank, L.R., Rand, T., Yeh, L., Wong, E.C., Skaf, A., Dantas, R.W., Haghighi, P., Trudell, D., Resnick, D., 1999. Collateral ligaments of the ankle: high-resolution MR imaging with a local gradient coil and anatomic correlation in cadavers 1. Radiographics 19, 673–683. Nie, B., Panzer, M.B., Mane, A., Mait, A., Donlon, J.P., Forman, J., Kent, R., 2016. A framework for parametric modeling of ankle ligaments to determine the in situ response under gross foot motion. Comput. Methods Biomech. Biomed. Eng. 19, 1254–1265. Nigg, B.M., Skarvan, G., Frank, C.B., Yeadon, M.R., 1990. Elongation and forces of ankle ligaments in a physiological range of motion. Foot Ankle Int. 11, 30–40. Ozeki, S., Yasuda, K., Kaneda, K., Yamakoshi, K., Yamanoi, T., 2002. Simultaneous strain measurement with determination of a zero strain reference for the medial and lateral ligaments of the ankle. Foot Ankle Int. 23, 825–832. Pie´rard, G.E., Lapie`re, C.M., 1987. Microanatomy of the dermis in relation to relaxed skin tension lines and Langer’s lines. Am. J. Dermatopathol. 9, 219–224. Rasmussen, O., 1985. Stability of the ankle joint. analysis of the function and traumatology of the ankle ligaments. Acta Orthop. Scand. 211, 1–75. Ray, R.G., Christensen, J.C., Gusman, D.N., 1997. Critical evaluation of anterior drawer measurement methods in the ankle. Clin. Orthop. Relat. Res. 334, 215–224. Reggiani, B., Leardini, A., Corazza, F., Taylor, M., 2006. Finite element analysis of a total ankle replacement during the stance phase of gait. J. Biomech. 39, 1435–1443. Roaas, A., Andersson, G.B., 1982. Normal range of motion of the hip, knee and ankle joints in male subjects, 30–40 years of age. Acta Orthop. 53, 205–208. Rubin, A., Sallis, R., 1996. Evaluation and diagnosis of ankle injuries. Am. Fam. Physician 54, 1609–1618. Shin, J., Yue, N., Untaroiu, C.D., 2012. A finite element model of the foot and ankle for automotive impact applications. Ann. Biomed. Eng. 40, 2519–2531. Siegler, S., Block, J., Schneck, C.D., 1988. The mechanical characteristics of the collateral ligaments of the human ankle joint. Foot Ankle Int. 8, 234–242. Simpson, T.W., Mauery, T.M., Korte, J.J., Mistree, F., 2001. Kriging models for global approximation in simulation-based multidisciplinary design optimization. AIAA J. 39, 2233–2241.

Subit, D., Masson, C., Brunet, C., Chabrand, P., 2008. Microstructure of the ligament-to-bone attachment complex in the human knee joint. J. Mech. Behav. Biomed. Mater. 1, 360–367. Takai, S., Woo, S.L.Y., Livesay, G.A., Adams, D.J., Fu, F.H., 1993. Determination of the in situ loads on the human anterior cruciate ligament. J. Orthop. Res. Off. Publ. Orthop. Res. Soc. 11, 686–695. Tannous, R.E., Bandok, F.A., Toridis, T.G., Eppinger, R.H., 1996. A three-dimensional finite element model of the human ankle: development and preliminary application to axial impulsive loading. In: Proceedings of the 40th Stapp Car Crash Conference. Society of Automotive Engineers (SAE), Warrendale, PA, USA, pp. 219–238. Viidik, A., 1990. Structure and Function of Normal and Healing Tendons and Ligaments, Biomechanics of Diarthrodial Joints. Springer-Verlag, New York3–38. Wei, F., Hunley, S.C., Powell, J.W., Haut, R.C., 2011. Development and validation of a computational model to study the effect of foot constraint on ankle injury due to external rotation. Ann. Biomed. Eng. 39, 756–765. Weindel, S., Schmidt, R., Rammelt, S., Claes, L., Campe, A.v., Rein, S., 2010. Subtalar instability: a biomechanical cadaver study. Arch. Orthop. trauma Surg. 130, 313–319. Weiss, J.A., Gardiner, J.C., 2001. Computational modeling of ligament mechanics. Crit. Rev. Biomed. Eng. 29, 1–70. Williams, G.N., Jones, M.H., Amendola, A., 2007. Syndesmotic ankle sprains in athletes. Am. J. Sports Med. 35, 1197–1207. Woo, S.L., Debski, R.E., Withrow, J.D., Janaushek, M.A., 1999. Biomechanics of knee ligaments. Am. J. Sports Med. 27, 533–543. Xenos, J.S., Hopkinson, W.J., Mulligan, M.E., Olson, E.J., Popovic, N. A., 1995. The tibiofibular syndesmosis. Evaluation of the ligamentous structures, methods of fixation, and radiographic assessment. J. Bone Jt. Surg. Am. 77, 847–856. Yoganandan, N., Kumaresan, S., Pintar, F.A., 2001. biomechanics of the cervical spine Part 2. Cervical spine soft tissue responses and biomechanical modeling. Clin. Biomech. 16, 1–27. Zhang, M., Meng, W., Davies, T.C., Zhang, Y., Xie, S.Q., 2016. A robot-driven computational model for estimating passive ankle torque with subject-specific adaptation. IEEE Trans. Biomed. Eng. 63, 814–821.