- 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)
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.
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].
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.
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).
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
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
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).
The predicted IDP was compared with the in vitro experimental data reported by Liu et al. [63], during 1.5 Nm flexion and extension
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).
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
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.
Torg JS, Guille JT, Jaffe S. Injuries to the cervical spine in American football players. J Bone Joint Surg Am. 2002;84:112–22.
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.
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.
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.
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.
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.
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.
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.
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.
Ö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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
Natarajan RN, Williams JR, Andersson GBJ. Recent advances in analytical modeling of lumbar disc degeneration. Spine (Phila Pa 1976). 2004;29:2733–41.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
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.
Reilly DT, Burstein AH, Frankel VH. The elastic modulus for bone. J Biomech. 1974;7(3):271.
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.
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.
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.
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.
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.
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.
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.
Yang KH, Kish VL. Compressibility measurement of human intervertebral nucleus pulposus. J Biomech Biomech. 1988;21(10):865.
Kleinberger M. Application of Finite Element Techniques to the Study of Cervical Spine Mechanics. In: 37th Stapp Car Crash Conference; 1993. p. 261–72.
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.
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.
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.
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.
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.
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.
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.
Leatt P, Reilly T, Troup JGD. Spinal loading during circuit weight-training and running. JSports Med. 1986;20(3):119–24.
Adams MA. Biomechanics of back pain. Acupunct Med Br Med Acupunct Soc. 2004;22:178–88.
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
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
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.
About this article
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
Received:
Accepted:
Published:
DOI: https://doi.org/10.1186/s12891-021-04172-1