Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

The importance of intervertebral disc material model on the prediction of mechanical function of the cervical spine

BMC Musculoskeletal Disorders volume 22, Article number: 324 (2021) Cite this article

Abstract

Background

Linear elastic, hyperelastic, and multiphasic material constitutive models are frequently used for spinal intervertebral disc simulations. While the characteristics of each model are known, their effect on spine mechanical response requires a careful investigation. The use of advanced material models may not be applicable when material constants are not available, model convergence is unlikely, and computational time is a concern. On the other hand, poor estimations of tissue’s mechanical response are likely if the spine model is oversimplified. In this study, discrepancies in load response introduced by material models will be investigated.

Methods

Three fiber-reinforced C2-C3 disc models were developed with linear elastic, hyperelastic, and biphasic behaviors. Three different loading modes were investigated: compression, flexion and extension in quasi-static and dynamic conditions. The deformed disc height, disc fluid pressure, range of motion, and stresses were compared.

Results

Results indicated that the intervertebral disc material model has a strong effect on load-sharing and disc height change when compression and flexion were applied. The predicted mechanical response of three models under extension had less discrepancy than its counterparts under flexion and compression. The fluid-solid interaction showed more relevance in dynamic than quasi-static loading conditions. The fiber-reinforced linear elastic and hyperelastic material models underestimated the load-sharing of the intervertebral disc annular collagen fibers.

Conclusion

This study confirmed the central role of the disc fluid pressure in spinal load-sharing and highlighted loading conditions where linear elastic and hyperelastic models predicted energy distribution different than that of the biphasic model.

Peer Review reports

Introduction

Among major spine injuries, cervical injuries could be life-altering as they might cause permanent loss of neural function. Cervical fractures are common in intensive recreational activities, such as American football and rugby [1, 2]. Experimental measurement provided invaluable insight into spine function to some extent, yet fails to quantify precisely other important parameters, such as strain/stress distribution under impact loads. Moreover, empirical investigations of some what-if questions would be time- and cost-expensive. Simulation modelling proved to be a viable means in biomechanical studies, in vivo or otherwise.

The spine simulation models have been improved extensively using calibrated constitutive models [3,4,5,6,7], detailed geometry reconstruction [8,9,10,11], and realistic physiological loading measurements [12,13,14]. While there is a general agreement on modelling the bony parts, various material models for the intervertebral disc (IVD) and ligaments were shown to be valid under certain loading conditions. Effects of ligament material modelling on spine mechanical behavior have been investigated extensively [8, 15,16,17,18]. However, studies on spine biomechanical response obtained from different IVD material models are sparse [19]. The linear elastic (LE), hyperelastic (HE), nucleus cavity and multiphasic theories are the most commonly used constitutive theories for simulating nucleus pulposus (NP) and annulus fibrosus (AF). The annular collagen fibers were simulated using spring elements with an organized spatial distribution in concentric lamellae [20]. The LE theory is the simplest material model to apply but is incapable of describing the nonlinear mechanical behaviour of IVD [3, 21]. The HE constitutive models addressed nonlinear stress-strain behaviour [22, 23], but the energy conservation nature of the LE and HE models failed to simulate the IVD creep/relaxation behaviour caused by fluid-solid interaction [24]. Studies that used LE and HE models estimated the intradiscal pressure (IDP) with the disc hydrostatic stress [25, 26]. The cavity model assumed the NP as a fluid volume enclosed within a membrane [27,28,29]. Although this assumption represented the load-bearing characteristic of the IDP under high loading rates, it ignored the shear stresses and failed to model progressive disc deformations under cyclic loads. The multiphasic theory [30, 31] was the most recent material model that showed the evident link between the fluid and solid interactions [19, 32, 33]. The viscoelastic model would also simulate a rate-dependent mechanical response [34], yet it lacks an accurate prediction of some tissue biomechanical phenomena, such as tissue swelling and degeneration.

The effects of LE, HE, and multiphasic IVD material models on mechanical response were studied using simple axisymmetric models [35], however, it was shown that spine components, i.e. ligaments, endplates, and facets, have a considerable contribution to the mechanical response of the functional spinal unit (FSU) [8, 32, 35, 36]. Nowadays, the spine components are reconstructed from computer tomography scans or magnetic resonance images with a high spatial resolution to predict more precise load-sharing of the FSU and simulate in-vivo loading conditions. Flexion, extension, lateral bending and axial compression are typical loading conditions supported by the cervical spine during daily activities [37,38,39,40,41]. Range of motion (RoM), facets contact pressure, disc height change, IDP, and creep/relaxation were the most important mechanical parameters of the FSU studied in the literature. Therefore, in the present study, a detailed nonlinear 3D finite element (FE) model of C2-C3 FSU was presented with three different fiber-reinforced IVD material models, including LE, HE, and Biphasic theories. The objective of the present study was to compare the kinetics and kinematics predicted by these models under flexion/extension moments and cyclic compression load.

Methodology

FE meshes and material models

A 3D FE model of the C2-C3 FSU was created (Fig. 1) and analyzed in Abaqus implicit software package (Dassault Systemes, 2017). Bony components geometry was reconstructed from 1 mm interval CT-Scans of a 39-year-old male [42]. These images were provided by the National Library of Medicine’s Visible Human Project (http://www.nlm.nih.gov/research/visible/visible_human.html). Each vertebra comprises cancellous and cortical bones, endplates, and articular facets (Fig. 1). The complex geometry of cortical bone was meshed using 3-node triangular shell elements of 1 mm thickness. Four-node tetrahedral solid elements were used to mesh the cancellous bone. Facet cartilages and endplates of 0.7 mm thickness were meshed with 8-node brick elements. Two capsules were created to connect the adjacent articular facets on both sides (Fig. 1) [22]. The IVD volume was divided into NP and AF ground substance with a proportion of 56 and 44%, respectively, according to the histological findings [43]. The AF mesh was generated by extruding seven layers of 8-node solid elements between two adjacent endplates, and reinforced in the circumferential direction by seven layers of collagen fibers (Fig. 1). Unidirectional nonlinear springs resisting tensile load only and organized in concentric lamellae with crosswise patterns close to ±35° were used to model the annular fibers [20]. The Anterior Longitudinal (ALL), Capsular (CL), Interspinous (ISL), Intertransverse (ITL), Ligamentum Flavum (LF), Posterior Longitudinal (PLL), and Supraspinous (SSL) ligaments were constructed using 4-node shell elements of 1 mm thickness by connecting their attachment points obtained from anatomical and histological studies to surrounding vertebrae [44,45,46].

Fig. 1
figure 1

The simulated C2-C3 FSU components in the present study

Full size image

Strain rate-dependent elasto-plastic Johnson-Cook material model was assigned to bony components [47, 48] (Table 1). The ligaments were defined by the Prony series viscoelastic model [47, 52, 53]. The LE material properties were assigned to the facet cartilages, and a frictionless surface-to-surface hard contact was defined between the two adjacent facets [9, 32] (Table 2). The effect of IVD material properties on the overall behaviour of FSU as well as the local spinal tissues was investigated by assigning LE, HE, and Biphasic material properties to the AF, NP and endplates (Table 3). All loads were applied to the centroid of C2, while the C3 bottom was completely fixed.

Table 1 Material properties of bony components for C2-C3 FSU model
Full size table
Table 2 Material properties of ligaments, IVD collagen fibers, and facet cartilages for C2-C3 FSU model. E = Young Modulus, ν = Poisson Ratio
Full size table
Table 3 Material properties of IVD components in the LE, HE, and Biphasic models of the C2-C3 FSU [35]
Full size table

Loading conditions

Cyclic compression

A 200 N compression ramp in 20 s was applied to the centroid of C2 following the Skrzypiec et al. [59] loading condition, and the C3 vertebral body was fully fixed. The IVD stress predicted by the developed FE models was compared with experimental measurements [59]. Furthermore, the performances of LE, HE, and Biphasic models were compared under a periodical compressive load with an amplitude of 100 N and frequency of 0.5 Hz in the form of F = 50 + 50 cos(π(t − 1)). This loading range was shown to cover the approximate 46 ± 7 N weight of human head [60, 61]. The cyclic compression was applied for 11,000 cycles, which is in line with the estimated number of steps a soldier takes during a 19 km walk to gain the Expert Infantryman Badge [62].

Sagittal bending

Quasi-static bending moments with magnitudes of 1 and 1.5 Nm were applied to the centroid of C2 in the anterior (+ve x-axis) and posterior (−ve x-axis), representing flexion and extension in the sagittal plane, respectively. The 1–1.5 Nm flexion/extension moments were judged to be sufficient to produce physiologic motions in the actual intact conditions [38, 63]. This loading condition was applied to validate the model with experimental results [38, 63]. The effect of the IVD material model on the FSU load-sharing was studied under a flexion/extension moment of 1.5 Nm in 30 s. Strain Energy stored in the FSU components was used to investigate their load-sharing.

Results

Cyclic compression

The IVD stress profile obtained from the Biphasic model, including the IDP and compressive stress, showed a good agreement with experimental measurements [59] (Fig. 2). The linear region of the IVD stress profile in the Biphasic model results was longer compared to experimental results [59]. The hydrostatic stresses of the LE and HE models exhibited a similar stress profile but with magnitudes smaller than the experimental values. The performance of LE, HE, and Biphasic models under the cyclic compression was compared in terms of facet joints contact pressure (Fig. 3), and disc height variation (Fig. 4). The contact pressure area in the left facet joint was larger than the right one in all models (Fig. 3). The LE and HE models produced the same contact pressure profile at every cycle, while in the biphasic model, the final contact pressure profile was developed over the load repetitions (Fig. 3). Figure 4 illustrated the disc height variation during the cyclic loading. No permanent height decrease was observed in the results of the LE and HE models. However, the disc height decreased exponentially from one cycle to the next cycle in the Biphasic model. An exponential curve was fitted to the mean of disc height variations in Fig. 4a. The decaying oscillation of the disc height was illustrated in Fig. 4b. A schematic diagram of the hysteresis loop in the Biphasic model was illustrated in the force-displacement diagram of Fig. 4c, where the dissipated energy increased with repeated loads. The strain energy ratios of major load-bearing soft tissues, including the IVD, endplates, and facets ratios were illustrated in Fig. 5. The strain energy could be a representation of the load-sharing distribution. In the LE and HE models, {IVD, endplates, facets} contributions to load-sharing were {72, 5, 23%} and {79, 8, 13%}, respectively, which remained uniform from cycle 1 to the last cycle. However, in the Biphasic model, the load-sharing varied through cycles; facets’ contribution increased from 1% in cycle 1 to 12% in cycle 250. On the other hand, endplates’ contribution decreased from 40 to 17% from cycle 1 to 250, respectively. After cycle 250, 500 s of loading, the load-sharing distribution had less than 0.1% between two consecutive cycles in the Biphasic model (Fig. 5).

Fig. 2
figure 2

The IVD stress profile along the sagittal midline diameter of the disc. The results were compared with the digitized curve from Skrzypiec et al. study [59]. The x-axis is the normalized length of the centerline, where x = 0 and 1 representing the posterior and anterior AF, respectively. In LE and HE models, the hydrostatic stress represented the IDP stress. S33: the compressive stress in the direction of the applied load

Full size image
Fig. 3
figure 3

Facet joints’ contact pressure distributions predicted by the LE, HE, and Biphasic models at the peak of cycles 1, 50, and 250. The LE and HE models resulted in an approximately identical contact pressure distribution. The contact pressure of the Biphasic model increased with load repetition

Full size image
Fig. 4
figure 4

a The disc height variation over 600 s of cyclic compression obtained from the LE, HE, and Biphasic models. The disc had a purely sinusoidal height change in the LE and HE models. However, the disc height variation in the Biphasic model had a decaying sinusoidal oscillation with the amplitude of δ as shown in b about an exponential curve with the equation of ∆ =  − 0.3 exp (0.004t) − 0.37. c Representative force-disc height variation in cycles 2, 3, 50, 100. The dissipated energy (highlighted area) increased with the number of load repetitions. ∆Ei is the dissipated energy at cycle i

Full size image
Fig. 5
figure 5

The predicted strain energy ratio stored in the IVD, endplates, and facets during the cyclic compression

Full size image

Sagittal bending

The FSU’s RoM predicted by all three models during 1 Nm and 1.5 Nm quasi-static flexion and extension were illustrated in Fig. 6. The predicted RoM in extension was almost identical in the three models. However, the RoM in flexion was higher in LE than the HE and Biphasic model. Fluid pressure in the NP predicted by the Biphasic model fell within 1 standard deviation (SD) of the IDP reported by Liu et al. [63] (Fig. 7). The LE and HE models did not consider the fluid phase, therefore the IDP was estimated by the hydrostatic pressure in these models and compared with the experimental IDP values (Fig. 7).

Fig. 6
figure 6

RoM (Degree) predicted by the Biphasic, LE and HE FEA models was compared with existing in vitro experimental data reported by Panjabi et al. [38], (1 Nm) and Liu et al. [63], (1.5 Nm) during flexion and extension

Full size image
Fig. 7
figure 7

The predicted IDP was compared with the in vitro experimental data reported by Liu et al. [63], during 1.5 Nm flexion and extension

Full size image

The normalized strain energy ratio of each component represented its contribution to the FSU load-sharing. Major load-bearing components (ligaments, IVD, endplates and facets) were considered to study the effect of the IVD material model on the FSU load-sharing in Figs. 8 and 9. A 1.5 Nm flexion/extension was applied linearly in 30 s and kept for 10,000 s for creep (Fig. 8). The strain energy ratios within the first 0.5° of RoM were excluded because initial facet contact pressures introduced noises in results. During flexion, {ligaments, IVD, endplates, facets} carried {54, 40, 6, 0%}, {50, 42, 8, 0%}, and {60, 26, 14, 0%} of the FSU’s total strain energy in the LE, HE, and Biphasic models, respectively (Fig. 8). While the components’ contribution to load-sharing was uniform in the LE and HE model, they redistributed considerably in the Biphasic model during the creep phase. Facets did not show any contribution during the flexion (Fig. 8).

Fig. 8
figure 8

Strain energy distribution among FSU components predicted by the LE, HE, and Biphasic material models during 1.5 Nm flexion/extension that was applied in 30 s (loading phase), and kept for 10,000 s (creep phase). The solid dot (●) shows the start time point of the creep phase in the Biphasic model

Full size image
Fig. 9
figure 9

The strain energy contribution of AF, NP and collagen fibers during 1.5 Nm flexion and extension, applied in 30 s and kept for 10,000 s. The solid dot (●) represents the start time point of the creep phase at t = 30 s

Full size image

During extension, the strain energy ratio in the ligaments and IVD decreased with further rotations, while load-sharing contributions of the endplates and facets increased (Fig. 8). This load-sharing redistribution behaviour was observed in all models under extension with almost identical magnitudes of {34, 25, 30, 11%}.

Within the IVD, the AF had the largest contribution to the load-bearing during flexion (in the range of 60–80%), following by the NP and collagen fibers (Fig. 9). Under flexion, the IVD was compressed almost everywhere in the LE and HE models, resulting in a negligible collagen fibers contribution to the load-sharing (less than 1%). However, in the Biphasic model, a group of collagen fibers were extended due to the IDP and carried about 15% of the IVD strain energy (Fig. 9. right). During extension, collagen fibers’ contribution to IVD mechanical response reached up to 52% with rotation. Similar to Fig. 8., a redistribution of load-sharing occurred in the IVD components during the creep phase for the Biphasic model. In flexion, the load-sharing was transferred from collagen fibers to the AF and NP, increasing their contribution from 60 to 75% and 17 to 26%, respectively, at the creep phase. Whereas no change was observed in the LE and HE results during the loading and creep phases (Fig. 9). The AF and NP had more contribution to flexion than extension in all models.

Discussion

The effect of LE, HE, and Biphasic IVD material models on the kinematics and kinetics responses of a C2-C3 model were investigated under flexion, extension, and cyclic compression. The solid phase in the LE and Biphasic models had identical material properties (Table 3), so any discrepancy between their mechanical responses would be caused by the solid-fluid interaction. The HE model was used to investigate whether the fluid-solid interaction effect in the FSU mechanical response could be described by material nonlinearity. The IVD material model affected the load-sharing, IDP and RoM under flexion and cyclic loadings more than the extension.

Before comparing the three LE, HE, and Biphasic models, their results were validated with experimental works for compression and flexion/extension loading conditions. The predicted IVD stresses under the 200 N ramp compression were compared with the experimental measurements in Fig. 2. The peak pressure at the NP and AF interface could be due to the change in permeability and elasticity modulus. From the linear part of the curve, the size of NP in [59] was smaller than the simulated NP in the present work. The disc pressure was approximately constant in the NP, and the mechanical response was governed mainly by the fluid pressure. This observation justified the use of cavity models in some studies to simulate the NP [27, 47, 64]. The predicted IDP in the NP under flexion and extension showed good agreement with the in vitro data measured by pressure transducers [63] (Fig. 7). The hydrostatic pressure of the NP in LE and HE models may not be a good representative of IDP in some loading conditions, such as extension (Fig. 7), however, it would be the most relevant parameter for IDP estimation in single-phase models. Peak IDP in extension was smaller than the flexion (Fig. 7). The FSU’ RoM predicted by all three models during 1 Nm and 1.5 Nm quasi-static flexion and extension were in good agreement with in vitro data [38, 63] except for the LE model during 1.5 Nm flexion, which was within two SD from the mean of experimental data (Fig. 6). The LE, HE, and Biphasic models predicted approximately identical RoMs (less than 2% difference) in extension (Fig. 6 - left), but their predictions diverged when the FSU was under flexion (Fig. 6 - right).

A cyclic compression was applied because it is a physiologically relevant loading condition. Moreover, under this dynamic loading condition, it was expected to obtain a clear discrepancy in the mechanical responses of the three IVD models. Overall, a repeated mechanical response was observed for the LE and HE models as these two models did not include any energy-dissipating element. Note that the viscoelastic ligaments were active only in extension. Contours of facets contact pressure (Fig. 3) during cyclic loads in the Biphasic model illustrated a growing load-bearing contribution from facet joints, which occurred mainly in the left facet. As a result, the contribution of facets in the strain energy increased from 1 to 12% in the Biphasic model after 250 cycles (Fig. 5). However, the LE and HE models resulted in an identical contact pressure contour at each cycle. The LE model was slightly softer than the HE model and resulted in about 12% more disc compression, as seen in Fig. 4a.

During a single compression cycle, the disc height variation was smaller in the Biphasic compared to the LE and HE models, indicating the importance of the IDP contribution to preventing excessive disc deformation under dynamic loading conditions. The LE/HE models predicted 0.32 and 0.27 mm disc heigh reduction at every cycle (Fig. 4a), respectively, while the disc heigh varied less than 0.1 mm at every cycle of compression (Fig. 4b). In the Biphasic model, the IVD did not have enough time during the unloading phase of cyclic compressions to recover its lost fluid, leading to an incremental disc height loss at every compression cycle that described a creep-like phenomenon. The exponential regression curve (R2 = 0.98) predicted a 0.37 mm permanent disc deformation after an infinite number of cyclic compressions. A corresponding value of 0.21 mm was obtained in an experimental study by Yoganandan et al. [65] (C4-C5 FSU, load amplitude = 150 N, f = 2 Hz). This permanent deformation in one spinal disc, when is cascaded to all vertebrae, may result in a 2–3 mm height loss [66]. The disc height variation had a decaying oscillation amplitude (δ) with time, started from δ = 0.09 mm in cycle 1 and decayed to a constant δ = 0.05 mm in cycle 250. The LE and HE models did not precisely predict the above discussed mechanical responses. Therefore, dynamic compression is an exemplary condition where the energy-conservative constitutive material models would fail to describe the IVD mechanical behaviour accurately.

Understanding spinal load-sharing is crucial in clinical studies [67]. In the present study, the strain energy of the FSU components was used to investigate the effect of IVD material models on the load-sharing response of FSU under flexion, extension and compression. In flexion, CL and PLL ligaments were the major load-bearing components in all three models, followed by IVD and endplates (Fig. 8). The ligaments and endplates strain energy ratios in the Biphasic model were about 10% greater than the corresponding values in LE and HE models. Reversely, the IVD strain energy ratio in the Biphasic model was 17% less than LE and HE models under flexion (Fig. 8). Facets contribution was negligible in flexion, yet reached 11% of total load-sharing in extension in all models (Fig. 8). The strain energy ratios were redistributed when creep started at t = 30 s in the Biphasic model (marked by a dot in Fig. 8); in extension, the role of facets in load-bearing increased, while the contribution of endplates and IVD decreased. This makes sense because when fluid exudes the IVD during the creep, it becomes deflated, and as a result, facets continue to carry more loads until their strain energy values reach the LE response at the end of the creep phase. The LE and HE models did not demonstrate this phenomenon.

The IVD components play distinct roles in various loading conditions. The AF was the major load-bearer in flexion, followed by NP and collagen fibers in all three models (Fig. 9). In extension, the AF was also the major load bearer at the beginning of loading, but with further extension, the collagen fibers showed more contribution than the AF. One major difference in the load-sharing response of the three models during flexion was the contribution of IVD collagen fibers; the LE and HE models predicted less than 2% strain energy ratio for collagen fibers, but this value was 20% in the Biphasic model. A higher contribution from collagen fibers under flexion reduced the contribution of AF from about 75% in the LE to 60% in the Biphasic model (Fig. 9. left). The IDP could cause the IVD volume increase in the Biphasic model and stretched collagen fibers during the loading phase. This could be considered an advantage of the Biphasic model, where an effective synergy among the IVD components was predicted which caused more uniform load-sharing distribution among components and reduced stresses. The discrepancy of load-sharing in the three models was less in extension than flexion. Apparently, in the creep phase, the IDP decreased gradually, and collagen fibers’ contribution reduced to about 2% as in the LE response. These results suggest that simulating a Biphasic IVD model is more important for flexion and compression than extension loading. A further investigation on considering the combined flexion/extension, lateral bending and axial compression may identify conditions that single-phase models would underperform. There are limitations to the presented FE simulation. Active musculature that exerts stabilizing forces on the spine was not simulated. The resultant force of spine muscles was summed by the applied loads in the present study. The series of multiple FSUs in the spine were simulated using an isolated C2-C3 FSU model, and the RoM was predicted within one SD from the mean of experimental measurements (Figs. 6 and 7). However, to simulate in vivo conditions more closely, a full spine model could be considered. Results of this study may be subject-specific; hence, sensitivity and parametric analyses could be conducted to investigate the reliability of the findings for different FSU sizes and material properties.

Conclusion

The single-phase IVD material models failed to describe the long-term IVD height decrease under cyclic loading. The contact area and contact pressure of facet joints were different in the Biphasic and LE/HE models. While a small difference was observed in the mechanical response of C2-C3 FSU under extension, the Biphasic and LE/HE models predicted distinct load-sharing during the loading phase of flexion. However, when the static condition is reached, Biphasic and LE models predicated the same mechanical response. These findings include that applying simpler material models, such as linear elastic, is justifiable when the spine is under quasi-static extension. However, more accurate material models, such as biphasic, are required for modeling intervertebral disc when the spine is under dynamic compression and flexion.

Availability of data and materials

Data sharing is not applicable to this article as no datasets were generated or analysed during the current study. The spine model used and analyzed in the current study is available upon reasonable request to Dr. Amin Komeili (AK; akomeili@uoguelph.ca).

Abbreviations

LE:

Linear Elastic

HE:

Hyperelastic

NP:

Nucleus Pulposus

AF:

Annulus Fibrosus

IVD:

Intervertebral Disc

IDP:

Intradiscal Pressure

FSU:

Functional Spinal Unit

RoM:

Range of Motion

FE:

Finite Element

CL:

Capsular Ligament

ISL:

Interspinous Ligament

ITL:

Intertransverse Ligament

LF:

Ligamentum Flavum Ligament

PLL:

Posterior Longitudinal Ligament

SSL:

Supraspinous Ligament

References

  1. Scher AT. Catastrophic rugby injuries of the spinal cord: changing patterns of injury. Br J Sports Med. 1991;25(1):57–60. https://doi.org/10.1136/bjsm.25.1.57.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  2. Torg JS, Guille JT, Jaffe S. Injuries to the cervical spine in American football players. J Bone Joint Surg Am. 2002;84:112–22.

    Article  Google Scholar 

  3. Holzapfel GA, Schulze-Bauer CAJ, Feigl G, Regitnig P. Single lamellar mechanics of the human lumbar anulus fibrosus. Biomech Model Mechanobiol. 2005;3(3):125–40. https://doi.org/10.1007/s10237-004-0053-8.

    Article  CAS  PubMed  Google Scholar 

  4. Schmidt H, Heuer F, Simon U, Kettler A, Rohlmann A, Claes L, Wilke HJ. Application of a new calibration method for a three-dimensional finite element model of a human lumbar annulus fibrosus. Clin Biomech. 2006;21(4):337–44. https://doi.org/10.1016/j.clinbiomech.2005.12.001.

    Article  Google Scholar 

  5. Yoganandan N, Kumaresan S, Pintar FA. Geometric and mechanical properties of human cervical spine ligaments. J Biomech Eng. 2000;122(6):623–9. https://doi.org/10.1115/1.1322034.

    Article  CAS  PubMed  Google Scholar 

  6. Shirazi-Adl A. On the fibre composite material models of disc annulus-comparison of predicted stresses. J Biomech. 1989;22(4):357–65. https://doi.org/10.1016/0021-9290(89)90050-X.

    Article  CAS  PubMed  Google Scholar 

  7. Kemper AR, McNally C, Duma SM. The influence of strain rate on the compressive stiffness properties of human lumbar intervertebral discs. Biomed Sci Instrum. 2007;43:176–81.

    PubMed  Google Scholar 

  8. Naserkhaki S, Arjmand N, Shirazi-Adl A, Farahmand F, El-Rich M. Effects of eight different ligament property datasets on biomechanics of a lumbar L4-L5 finite element model. J Biomech. 2018;70:33–42. https://doi.org/10.1016/j.jbiomech.2017.05.003.

    Article  CAS  PubMed  Google Scholar 

  9. Dreischarf M, Zander T, Shirazi-Adl A, Puttlitz CM, Adam CJ, Chen CS, Goel VK, Kiapour A, Kim YH, Labus KM, Little JP, Park WM, Wang YH, Wilke HJ, Rohlmann A, Schmidt H. Comparison of eight published static finite element models of the intact lumbar spine: predictive power of models improves when combined together. J Biomech. 2014;47(8):1757–66. https://doi.org/10.1016/j.jbiomech.2014.04.002.

    Article  CAS  PubMed  Google Scholar 

  10. Schmidt H, Heuer F, Drumm J, Klezl Z, Claes L, Wilke HJ. Application of a calibration method provides more realistic results for a finite element model of a lumbar spinal segment. Clin Biomech. 2007;22(4):377–84. https://doi.org/10.1016/j.clinbiomech.2006.11.008.

    Article  Google Scholar 

  11. Östh J, Brolin K, Svensson MY, Linder A. A Female Ligamentous Cervical Spine Finite Element Model Validated for Physiological Loads. J Biomech Eng. 2016;138(6):061005.

    Article  Google Scholar 

  12. Schmidt H, Shirazi-Adl A, Galbusera F, Wilke HJ. Response analysis of the lumbar spine during regular daily activities-a finite element analysis. J Biomech. 2010;43(10):1849–56. https://doi.org/10.1016/j.jbiomech.2010.03.035.

    Article  PubMed  Google Scholar 

  13. Rohlmann A, Graichen F, Kayser R, Bender A, Bergmann G. Loads on a telemeterized vertebral body replacement measured in two patients. Spine (Phila Pa 1976). 2008;33(11):1170–9. https://doi.org/10.1097/BRS.0b013e3181722d52.

    Article  Google Scholar 

  14. Wilke H, Neef P, Caimi M, Hoogland T, Claes LE. New In Vivo measurements of pressures in the intervertebral disc in daily life. Spine (Phila Pa 1976). 1999;24(8):755.

    Article  CAS  Google Scholar 

  15. Brolin K, Halldin P. Development of a finite element model of the upper cervical spine and a parameter study of ligament characteristics. Spine (Phila Pa 1976). 2004;29(4):376–85. https://doi.org/10.1097/01.BRS.0000090820.99182.2D.

    Article  Google Scholar 

  16. Mattucci SFE, Cronin DS. 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. 2015;41:251–60. https://doi.org/10.1016/j.jmbbm.2014.09.023.

    Article  PubMed  Google Scholar 

  17. Mattucci SFE, Moulton JA, Chandrashekar N, Cronin DS. Strain rate dependent properties of younger human cervical spine ligaments. J Mech Behav Biomed Mater. 2012;10:216–26. https://doi.org/10.1016/j.jmbbm.2012.02.004.

    Article  PubMed  Google Scholar 

  18. Sterba M, Aubin CÉ, Wagnac E, Fradet L, Arnoux PJ. Effect of impact velocity and ligament mechanical properties on lumbar spine injuries in posterior-anterior impact loading conditions: a finite element study. Med Biol Eng Comput. 2019;57(6):1381–92. https://doi.org/10.1007/s11517-019-01964-5.

    Article  PubMed  Google Scholar 

  19. Galbusera F, Schmidt H, Noailly J, Malandrino A, Lacroix D, Wilke HJ, Shirazi-Adl A. Comparison of four methods to simulate swelling in poroelastic finite element models of intervertebral discs. J Mech Behav Biomed Mater. 2011;4(7):1234–41. https://doi.org/10.1016/j.jmbbm.2011.04.008.

    Article  PubMed  Google Scholar 

  20. Naserkhaki S, Jaremko JL, Adeeb S, El-Rich M. On the load-sharing along the ligamentous lumbosacral spine in flexed and extended postures: finite element study. J Biomech. 2016;49(6):974–82. https://doi.org/10.1016/j.jbiomech.2015.09.050.

    Article  PubMed  Google Scholar 

  21. Iatridis JC, Setton LA, Foster RJ, Rawlins BA, Weidenbaum M, Mow VC. Degeneration affects the anisotropic and nonlinear behaviors of human anulus fibrosus in compression. J Biomech. 1998;31(6):535–44. https://doi.org/10.1016/S0021-9290(98)00046-3.

    Article  CAS  PubMed  Google Scholar 

  22. Mustafy T, Moglo K, Adeeb S, El-Rich M. Injury mechanisms of the ligamentous cervical C2-C3 functional spinal unit to complex loading modes: finite element study. J Mech Behav Biomed Mater. 2016;53:384–96. https://doi.org/10.1016/j.jmbbm.2015.08.042.

    Article  PubMed  Google Scholar 

  23. Arjmand N, Shirazi-Adl A, Parnianpour M. Trunk biomechanical models based on equilibrium at a single-level violate equilibrium at other levels. Eur Spine J. 2007;16(5):701–9. https://doi.org/10.1007/s00586-006-0263-0.

    Article  CAS  PubMed  Google Scholar 

  24. Komeili A, Abusara Z, Federico S, Herzog W. Effect of strain rate on transient local strain variations in articular cartilage. J Mech Behav Biomed Mater. 2019;95:60–6. https://doi.org/10.1016/j.jmbbm.2019.03.022.

    Article  PubMed  Google Scholar 

  25. Brolin K, Hedenstierna S, Halldin P, Bass C, Alem N. The importance of muscle tension on the outcome of impacts with a major vertical component. Int J Crashworthiness. 2008;13(5):487–98. https://doi.org/10.1080/13588260802215510.

    Article  Google Scholar 

  26. Cai XY, YuChi CX, Du CF, Mo ZJ. The effect of follower load on the range of motion, facet joint force, and intradiscal pressure of the cervical spine: a finite element study. Med Biol Eng Comput. 2020;58(8):1695–705. https://doi.org/10.1007/s11517-020-02189-7.

    Article  PubMed  Google Scholar 

  27. Khoddam-Khorasani P, Arjmand N, Shirazi-Adl A. Effect of changes in the lumbar posture in lifting on trunk muscle and spinal loads: a combined in vivo, musculoskeletal, and finite element model study. J Biomech. 2020;104:109728. https://doi.org/10.1016/j.jbiomech.2020.109728.

    Article  CAS  PubMed  Google Scholar 

  28. Eskandari AH, Arjmand N, Shirazi-Adl A, Farahmand F. Hypersensitivity of trunk biomechanical model predictions to errors in image-based kinematics when using fully displacement-control techniques. J Biomech. 2019;84:161–71. https://doi.org/10.1016/j.jbiomech.2018.12.043.

    Article  CAS  PubMed  Google Scholar 

  29. Park WM, Kim K, Kim YH. Effects of degenerated intervertebral discs on intersegmental rotations, intradiscal pressures, and facet joint forces of the whole lumbar spine. Comput Biol Med. 2013;43(9):1234–40. https://doi.org/10.1016/j.compbiomed.2013.06.011.

    Article  PubMed  Google Scholar 

  30. Ehlers W, Karajan N, Markert B. A porous media model describing the inhomogeneous behaviour of the human intervertebral disc. Materwiss Werksttech. 2006;37(6):546–51. https://doi.org/10.1002/mawe.200600035.

    Article  CAS  Google Scholar 

  31. Ehlers W, Karajan N, Markert B. An extended biphasic model for charged hydrated tissues with application to the intervertebral disc. Biomech Model Mechanobiol. 2009;8(3):233–51. https://doi.org/10.1007/s10237-008-0129-y.

    Article  CAS  PubMed  Google Scholar 

  32. Williams JR, Natarajan RN, Andersson GBJ. Inclusion of regional poroelastic material properties better predicts biomechanical behavior of lumbar discs subjected to dynamic loading. J Biomech. 2007;40(9):1981–7. https://doi.org/10.1016/j.jbiomech.2006.09.022.

    Article  PubMed  Google Scholar 

  33. Natarajan RN, Williams JR, Andersson GBJ. Recent advances in analytical modeling of lumbar disc degeneration. Spine (Phila Pa 1976). 2004;29:2733–41.

    Article  Google Scholar 

  34. Barker JB, Cronin DS, Nightingale RW. Lower cervical spine motion segment computational model validation: kinematic and kinetic response for quasistatic and dynamic loading. J Biomech Eng. 2017;139:6.

    Article  Google Scholar 

  35. Nikkhoo M, Haghpanahi M, Parnianpour M, Wang JL. An axisymmetric poroelastic model for description of the short-term and long-term creep behavior of L4-L5 intervertebral disc. In: 2011 1st Middle East Conference on Biomedical Engineering, MECBME 2011; 2011. p. 308–11.

    Google Scholar 

  36. Malandrino A, Planell JA, Lacroix D. Statistical factorial analysis on the poroelastic material properties sensitivity of the lumbar intervertebral disc under compression, flexion and axial rotation. J Biomech. 2009;42(16):2780–8. https://doi.org/10.1016/j.jbiomech.2009.07.039.

    Article  PubMed  Google Scholar 

  37. Panjabi MM, Summers DJ, Pelker RR, Videman T, Friedlaender GE, Southwick WO. Three-dimensional load-displacement curves due to froces on the cervical spine. J Orthop Res. 1986;4(2):152–61. https://doi.org/10.1002/jor.1100040203.

    Article  CAS  PubMed  Google Scholar 

  38. Panjabi MM, Crisco JJ, Vasavada A, Oda T, Cholewicki J, Nibu K, Shin E. Mechanical properties of the human cervical spine as shown by three-dimensional load-displacement curves. Spine (Phila Pa 1976). 2001;26(24):2692–700. https://doi.org/10.1097/00007632-200112150-00012.

    Article  CAS  Google Scholar 

  39. Bell KM, Yan Y, Debski RE, Sowa GA, Kang JD, Tashman S. Influence of varying compressive loading methods on physiologic motion patterns in the cervical spine. J Biomech. 2016;49(2):167–72. https://doi.org/10.1016/j.jbiomech.2015.11.045.

    Article  PubMed  Google Scholar 

  40. Emanuel KS, van der Veen AJ, Rustenburg CME, Smit TH, Kingma I. Osmosis and viscoelasticity both contribute to time-dependent behaviour of the intervertebral disc under compressive load: a caprine in vitro study. J Biomech. 2018;70:10–5. https://doi.org/10.1016/j.jbiomech.2017.10.010.

    Article  PubMed  Google Scholar 

  41. Bell K, Yan Y, Dugan T, Kang J. Effect of follower load application on intact cervical spine kinetics. In: ORS 2013 Annual Meeting. Pittsburgh; 2013.

    Google Scholar 

  42. Moghaddam F. A 3D continuum finite element muscle model for the investigation of cervical spine load-sharing mechanisms and injury assessment during impact loading scenarios by Fatemeh Moghaddam; 2018.

    Google Scholar 

  43. Panjabi MM, Oxland TR, Parks EH. Quantitative anatomy of cervical spine ligaments. Part ii. Middle and lower cervical spine. J Spinal Disord. 1991;4(3):278–85.

    Google Scholar 

  44. Yoganandan N, Kumaresan S, Pintar FA. Biomechanics of the cervical spine. Part 2. Cervical spine soft tissue responses and biomechanical modeling. Clin Biomech Elsevier Science Ltd. 2001;16:1–27.

    Article  CAS  Google Scholar 

  45. Panjabi M, Dvorak J, Crisco JJ, Oda T, Hilibrand A, Grob D. Flexion, extension, and lateral bending of the upper cervical spine in response to alar ligament transections. J Spinal Disord. 1991;4(2):157–67. https://doi.org/10.1097/00002517-199106000-00005.

    Article  CAS  PubMed  Google Scholar 

  46. Panjabi M, Dvorak J, Crisco JJ, Oda T, Wang P, Grob D. Effects of alar ligament transection on upper cervical spine rotation. J Orthop Res. 1991;9(4):584–93. https://doi.org/10.1002/jor.1100090415.

    Article  CAS  PubMed  Google Scholar 

  47. Mustafy T, El-Rich M, Mesfar W, Moglo K. Investigation of impact loading rate effects on the ligamentous cervical spinal load-partitioning using finite element model of functional spinal unit C2-C3. J Biomech. 2014;47(12):2891–903. https://doi.org/10.1016/j.jbiomech.2014.07.016.

    Article  PubMed  Google Scholar 

  48. El-Rich M, Arnoux PJ, Wagnac E, Brunet C, Aubin CE. Finite element investigation of the loading rate effect on the spinal load-sharing changes under impact conditions. J Biomech. 2009;42(9):1252–62. https://doi.org/10.1016/j.jbiomech.2009.03.036.

    Article  PubMed  Google Scholar 

  49. Reilly DT, Burstein AH, Frankel VH. The elastic modulus for bone. J Biomech. 1974;7(3):271.

    Article  CAS  Google Scholar 

  50. DeWit JA, Cronin DS. Cervical spine segment finite element model for traumatic injury prediction. J Mech Behav Biomed Mater. 2012;10:138–50. https://doi.org/10.1016/j.jmbbm.2012.02.015.

    Article  PubMed  Google Scholar 

  51. Kopperdahl DL, Keaveny TM. Yield strain behavior of trabecular bone. J Biomech. 1998;31(7):601–8. https://doi.org/10.1016/S0021-9290(98)00057-8.

    Article  CAS  PubMed  Google Scholar 

  52. Wagnac E, Arnoux P-J, Anaïs G, Aubin C-E. Finite element analysis of the influence of loading rate on a model of the full lumbar spine under dynamic loading conditions. Med Biol Eng Comput. 2012;50(9):903–15. https://doi.org/10.1007/s11517-012-0908-6.

    Article  PubMed  Google Scholar 

  53. Tchako A, Sadegh AM. Stress changes in intervertebral discs of the cervical spine after partial fusion. In: American Society of Mechanical Engineers, Bioengineering Division (Publication) BED. American Society of Mechanical Engineers (ASME); 2003. p. 375–6.

    Google Scholar 

  54. Sadegh AM, Abraham T. Vertebral stress of a cervical spine model under dynamic load. Technol Health Care. 2000;8(2):143–54. https://doi.org/10.3233/THC-2000-8205.

    Article  CAS  PubMed  Google Scholar 

  55. Wagnac E, Arnoux PJP-J, Garo A, Aubin C-EC, Anaïs G, et al. Finite element analysis of the influence of loading rate on a model of the full lumbar spine under dynamic loading conditions. Med Biol Eng Comput. 2012;50(9):903–15. https://doi.org/10.1007/s11517-012-0908-6.

    Article  PubMed  Google Scholar 

  56. Yoganandan N, Pintar FA, Maiman DJ, Cusick JF, Sances A, Walsh PR. Human head- neck biomechanics under axial tension. Med Eng Phys. 1996;18(4):289–94. https://doi.org/10.1016/1350-4533(95)00054-2.

    Article  CAS  PubMed  Google Scholar 

  57. Yang KH, Kish VL. Compressibility measurement of human intervertebral nucleus pulposus. J Biomech Biomech. 1988;21(10):865.

    Article  Google Scholar 

  58. Kleinberger M. Application of Finite Element Techniques to the Study of Cervical Spine Mechanics. In: 37th Stapp Car Crash Conference; 1993. p. 261–72.

    Google Scholar 

  59. Skrzypiec DM, Pollintine P, Przybyla A, Dolan P, Adams MA. The internal mechanical properties of cervical intervertebral discs as revealed by stress profilometry. Eur Spine J. 2007;16(10):1701–9. https://doi.org/10.1007/s00586-007-0458-z.

    Article  PubMed  PubMed Central  Google Scholar 

  60. Yoganandan N, Pintar FA, Zhang J, Baisden JL. Physical properties of the human head: Mass, center of gravity and moment of inertia. J Biomech Elsevier. 2009;42:1177–92.

    Article  Google Scholar 

  61. Plagenhoef S, Gaynor Evans F, Abdelnour T. Anatomical data for analyzing human motion. Res Q Exerc Sport. 1983;54(2):169–78. https://doi.org/10.1080/02701367.1983.10605290.

    Article  Google Scholar 

  62. Motiwale S, Subramani AV, Zhou X, Kraft RH. Damage prediction for a cervical spine intervertebral disc. In: ASME International Mechanical Engineering Congress and Exposition, Proceedings (IMECE). American Society of Mechanical Engineers (ASME); 2016.

    Book  Google Scholar 

  63. Liu Q, Guo Q, Yang J, Zhang P, Xu T, Cheng X, Chen J, Guan H, Ni B. Subaxial cervical Intradiscal pressure and segmental kinematics following Atlantoaxial fixation in different angles. World Neurosurg. 2016;87:521–8. https://doi.org/10.1016/j.wneu.2015.09.025.

    Article  PubMed  Google Scholar 

  64. Dehghan-Hamani I, Arjmand N, Shirazi-Adl A. Subject-specific loads on the lumbar spine in detailed finite element models scaled geometrically and kinematic-driven by radiography images. Int J Numer Method Biomed Eng. 2019;35(4):e3182. https://doi.org/10.1002/cnm.3182.

    Article  PubMed  Google Scholar 

  65. Yoganandan N, Umale S, Stemper B, Snyder B. Fatigue responses of the human cervical spine intervertebral discs. J Mech Behav Biomed Mater. 2017;69:30–8. https://doi.org/10.1016/j.jmbbm.2016.11.026.

    Article  PubMed  Google Scholar 

  66. Leatt P, Reilly T, Troup JGD. Spinal loading during circuit weight-training and running. JSports Med. 1986;20(3):119–24.

    CAS  Google Scholar 

  67. Adams MA. Biomechanics of back pain. Acupunct Med Br Med Acupunct Soc. 2004;22:178–88.

    Article  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

This work was supported by The Natural Sciences and Engineering Research Council of Canada (NSERC) [grant number 401611, 2020].

Author information

Authors and Affiliations

  1. School of Engineering, University of Guelph, Guelph, Canada

    Amin Komeili, Akbar Rasoulian & Fatemeh Moghaddam

  2. Healthcare Engineering Innovation Center, Department of Mechanical Engineering, Khalifa University, Abu Dhabi, United Arab Emirates

    Marwan El-Rich

  3. Department of Mechanical and Manufacturing Engineering, University of Calgary, Calgary, Canada

    Le Ping Li

Authors
  1. Amin Komeili
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Akbar Rasoulian
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Fatemeh Moghaddam
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Marwan El-Rich
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Le Ping Li
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

A shared co-first authorship is given to Drs. Komeili, Rasoulian and Moghaddam who equally contributed to this work. A.K conducted finite element analyses, extracted results, and wrote the main manuscript text. A.R helped in writing the introduction section, conducted the literature review, and prepared figures and tables. F.M developed the finite element model and meshed parts. M.E provided the spine model and revised the manuscript for important intellectual content. L.L provided insight into the technical matters and revised the manuscript. All authors read and approved the final version of the manuscript.

Corresponding author

Correspondence to Amin Komeili.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

Amin Komeili has no competing interests.

Akbar Rasoulian has no competing interests.

Fatemeh Moghaddam has no competing interests.

Marwan El-Rich has no competing interests.

LePing Li is a senior editorial board member of this journal.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Komeili, A., Rasoulian, A., Moghaddam, F. et al. The importance of intervertebral disc material model on the prediction of mechanical function of the cervical spine. BMC Musculoskelet Disord 22, 324 (2021). https://doi.org/10.1186/s12891-021-04172-1

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12891-021-04172-1

BMC Musculoskeletal Disorders

ISSN: 1471-2474

Contact us