The biomechanical study of cervical spine: A Finite Element Analysis

Abstract

The biomechanical study helps us to understand the mechanics of the human cervical spine. A three dimensional Finite Element (FE) model for C3 to C6 level was developed using computed tomography (CT) scan data to study the mechanical behaviour of the cervical spine. A moment of 1 Nm was applied at the top of C3 vertebral end plate and all degrees of freedom of bottom end plate of C6 were constrained. The physiological motion of the cervical spine was validated using published experimental and FE analysis results. The von Mises stress distribution across the intervertebral disc was calculated along with range of motion. It was observed that the predicted results of functional spine units using FE analysis replicate the real behaviour of the cervical spine.

Keywords

Cervical spine range of motion disc stress Finite Element Analysis biomechanical study patient specific

Introduction

The cervical spine is a complex structure in the anatomy of the human body. The intervertebral disc (IVD) consists of two main components: annulus fibrosus (AF) and nucleus pulposus (NP). The AF forms the outer part, it withstands large and complex loads. The AF is the major load-bearing component of the IVD.^1,2 The main function of NP is to dissipate the compressive load and acts like a shock absorber.³ A healthy disc is responsible for the flexibility and mobility of every individual spinal segment.⁴ The higher level of risk factor of the neck pain is the person becoming disabled.^5–7 The mechanical principles and the mechanics involved in the structure help the physicians to understand and plan the treatment strategies.^8–10 In vitro model, in vivo model and mathematical models such as finite element (FE) models are used to interpret the underlying mechanism of the cervical spine.¹¹ To study the underlying biomechanics, in vitro model has two existing techniques namely, animal cadaver testing and human cadaver testing.^12–14 The major drawbacks commonly observed in the study of human cervical spine are demographics based, anatomical based, structural based and physical properties variability. The in vivo study helps us to understand the kinematics and the muscle forces generated in the cervical spine.^15,16 The range of cervical disc degeneration analysed during dynamic neck motion as in Yu et al.¹⁷ could be used for understanding the cervical spine disc degeneration mechanics. An anatomically accurate FE model of the human cervical spine offers a better understanding in biomechanics and it needs to be validated with the experimental studies.^18–21 Also, the FE models helps to find the intervertebral disc mechanics, injury mechanism and contributes in the diagnosis of cervical spine problems and treatment methodologies.^22–25

The objective of this study is to predict the intervertebral stress values across C3–C4, C4–C5 and C5–C6 level for the flexion, extension, lateral bending right, axial rotation right using FEA.

Materials and methods

To create a three dimensional (3D) FE model^26,27 of cervical spine computed tomography (CT) scan of 35-year healthy male subject was used (Figure 1) The geometric data (66 slices, 0.5 mm slice thickness) of C3–C6 level was obtained from the CT data. The CT scan images were processed in an image processing software (Mimics 10.0; Materialise Technologies, Leuven, Belgium). By using segmentation, the surface model of the vertebral levels from C3 to C6 level were obtained. The model was meshed using the commercial FE package ANSYS 18.2 (ANSYS, Inc, USA) as shown in Figure 1.

Figure 1.

Intact FE model C3–C6.

During modelling of the vertebrae, the thickness of the cortical shell was kept as 0.5 mm.²⁸ The anterior sides of the vertebrae were meshed using eight noded brick elements and the posterior sides of the vertebra were meshed as tetrahedral elements due to the anatomical irregularity. The element type and the material properties were assigned as listed in Table 1.

Table 1.

Material properties of cervical spine model.

Component	Element type	Young’s modulus (MPa)	Poisson’s ratio	Cross-sectional area (mm²)	References
Bony Structures
Cortical bone		12,000	0.29		Zhang et al.²²
Cancellous bone	Solid	450	0.29		Zhang et al.²²
Posterior bone		3500	0.29		Zhang et al.²²
End plate		500	0.40		Zhang et al.²²
Intervertebral disc
Annulus fibrosus	Solid	3.4	0.40		Zhang et al.²²
Nucleus pulposus		1.0	0.49		Zhang et al.²²
Ligaments
ALL		30.0		6	Zhang et al.^20,22
PLL		20.0		5	Zhang et al.^20,22
ISL	Link (tension only)	1.5		10	Zhang et al.^20,22
LF		1.5		10	Zhang et al.^20,22
CL		20.0		5	Zhang et al.^20,22

The anterior side and the posterior size disc heights were maintained in the FE model.²⁹ The intervertebral discs were modelled as two separate components that is, annulus fibrosis and nucleus pulposus. The volume occupied by the annulus was around 60% and the nucleus was around 40%.³⁰ Endplates were modelled at the inferior and the superior sides of the vertebral bodies as 0.6 mm.³¹

The major ligaments in the cervical spine that is, ALL, anterior longitudinal ligament; PLL, posterior longitudinal ligament; ISL, interspinous ligament; LF, ligamentum flavum; CL, capsular ligament were modelled using two noded link elements.³² In this study, in order to replicate the real phenomenon of the ligament tension only link elements were considered. The material property and the cross sectional area of the ligaments were assigned as listed in Table 1.

The intervertebral disc was considered as the main load bearing area of the functional spine unit so the anterior side was meshed using hexahedral elements. The average mesh element size for FE model is 0.5 mm. Total elements of the model were 819,060 and convergence test done as per literature.³³

In the FE analysis, four loading conditions that is, flexion, extension, right lateral bending and right axial rotation were considered as shown in Figure 2. For each loading condition, 1 Nm moment was applied.³⁴ To create the moment, force was applied accordingly at the top surface of the C3 vertebra by creating a pilot node. All degree of freedom of bottom surface of the C6 level were completely constrained.

Figure 2.

Boundary condition.

Result

Validation

The C3–C6 FE model was analysed under static loading condition and the predicted ROM were compared with the results from the experimental studies^{21,23,36–39} and finite element analysis.^21,35,40 The physiological motions were validated with the experimental studies, based on this validation the von Mises stress distribution across the intervertebral discs were calculated.

Range of motion and the maximum von Mises stress

In the flexion condition, analysis was conducted for 1 Nm and the intersegmental motion across C3–C6 were calculated. During flexion, the ROM observed in C3–C4, C4–C5 and C5–C6 levels were 5.84°, 4.2° and 5.19° respectively. In extension, ROM found in C3–C4, C4–C5 and C5–C6 levels were 4.13°, 2.5° and 4.3° respectively. In lateral bending, the ROM observed in C3–C4, C4–C5 and C5–C6 levels were 3.1°, 3.4° and 2.8° respectively. In axial rotation, the ROM observed in C3–C4, C4–C5 and C5–C6 levels were 2.6°, 3.03° and 2.9° respectively.

The maximum von Mises stress was calculated for the C3–C4, C4–C5 and C5–C6 level discs. The von Mises stress distribution in the disc are shown in Figure 3. During flexion, the maximum von Mises stress across the C3–C4 disc was calculated as 0.57 MPa, in C4–C5 disc it was 0.51 MPa and in C5–C6 disc it was 0.51 MPa. During extension, the maximum stress observed in C3–C4, C4–C5 and C5–C6 discs were 0.57 MPa, 0.45 MPa and 0.51 MPa respectively. In comparison with the published results³⁴ the flexion exhibits more intervertebral disc stress in comparison with the extension. In lower flexion, the load carrying proportion of the disc is high. In case of higher flexion, the ligaments will play a major role.¹⁹ While in extension, the facets play a major role,¹⁹ thus the disc stress was less in extension compared to flexion. In case of lateral bending, the maximum stress found in the discs were 0.97 MPa in C3–C4, 0.47 MPa in C4–C5, and 0.65 MPa in C5–C6. Under axial rotation, observed maximum stresses were 2.04 MPa, 1.9 MPa and 1.6 MPa in C3–C4, C4–C5 and C5–C6 discs respectively. In lateral bending and axial rotation, the disc is the major load bearing area so the intervertebral disc stresses will be more in axial rotation than lateral bending.

Figure 3.

Summary of intervertebral disc stresses for C3–C6 for various physiological motion.

Discussion

Range of motion

The three dimensional FE model was created to study the physiological motion of the cervical spine from C3 to C6 level. The main advantage of this FE study is we compared our FE model with both in vitro and in silico study. In our study we used percentage to signify the range of motion when compared with the in silico studies the higher value indicates it is flexible and the lower value indicates it is stiffer. In addition to it, we verified our FE study range of motion with minimum value and maximum value of the in vitro studies.^23,36–39

In flexion as shown in the Figure 4(a) while we are comparing with the Cai et al.³⁵ our FE model behaves 12% higher and flexible in the C3–C4 level. In the C4–C5 level we found it is 5% less and a stiffer behaviour is seen. In the C5–C6 level we found it is 15% higher than Cai et al.³⁵ and a flexible behaviour is seen. In comparison with in vitro study Wheeldon et al.³⁷ in the C4–C5 level of our FE model exhibits stiffer values. In the flexion in vitro studies the motion of the C3–C4 level it is varying from 1.4° to 10.07 ° in C4–C5 level it varies from 2.29° to 8.54° and in C5–C6 level it is around 2.46° to 9.89°. In comparison with the published in vitro studies^23,36,37 our FE model behaves within the limits in all the segments.

Figure 4.

ROM predicted by intact FE model compared with experimental and FE studies (a) Flexion (b) Extension (c) Lateral bending (d) Axial rotation.

In the extension as shown in the Figure 4(b) the C3–C4 level in comparison with the Cai et al.³⁵ we found it is 29% higher and the FE model is flexible. With the C4–C5 level we found it is 14% less and stiffer and in the C5–C6 level it’s seen that its 6% higher than the published in silico studies.³⁵ In comparison with the in vitro studies it is observed that a good agreement is seen in the C3–C4 with^23,37 but it is flexible in comparison with Traynelis et al.³⁶ In the C4–C5 level it is stiffer in comparison with Panjabi et al.²³ and is within the limits with the published in vitro study.^36,37 In the C5–C6 level the values are within the limits as in Panjabi et al.²³ and Wheeldon et al.³⁷ and a flexibility is seen in comparison with Traynelis et al.³⁶ In the extension in vitro studies the motion of the C3–C4 level is varying from 1.32° to 5.5° in the C4–C5 level it varies from 0.75° to 6.71°and from C5–C6 level it is around 1.07° to 7.21° while comparing with in vitro studies^23,36,37 all the segments were within the limits.

In the lateral bending as shown in the Figure 4(c) in comparison with the in silico study Cai et al.³⁵ we found in the C3–C4 it is 28% in C4–C5 it is 70% and in C5–C6 it is 47% higher. In this case our FE model exhibits flexibility. With in vitro study we found in the C3–C4 level it is stiffer²³ while with Yoganandan et al.³⁸ it is within the limits. In the C4–C5 level it is stiffer compared with Panjabi et al.²³ whereas made a good agreement with Traynelis et al.³⁶ and Yoganandan et al.³⁸ In the C5–C6 level a good agreement is seen with all the in vitro studies.^23,36,38 In the lateral bending in vitro studies^23,36,38 in the C3–C4 level we found the value of bending is around 2.75° to 6.35° in the C4–C5 level it is around 0.6° to5.5° in the C5–C6 level we found it is around 1.6° to 5.28° our FE model value were within the limits.

In the axial rotation as shown in the Figure 4(d) in comparison with the Cai et al.³⁵ we found in all the levels our FE model exhibits a stiffer behaviour. In the C3–C4 level it is 12% lesser and in C4–C5 it 2% less and in the C5–C6 it is 4% less values were seen. In the axial rotation in all the segment we found a good agreement with the published in vitro studies.^23,36,39 In the axial rotation in C3–C4 level the value of the rotation is around 1.5° to 4.0° in the C4–C5 level it is around 1.60° to 3.92° and in the C5–C6 level it is around 1.85° to 3.96°. We found the FE model rotation value were within the limits. With Ha et al.⁴⁰ we found our FE model is flexible in all the physiological motion. This may be due the bone quality, material property assigned and they applied a preload with moment whereas we applied only pure moments. This may be the reason that our model is flexible. Our FE model in comparison with Kallemeyn et al.²¹ both in vitro and in silico study our FE model exhibits stiffer behaviour in all the physiological motion. This behaviour may be due to the reason they calibrated their FE model with the experimental results.

Stress distribution in the intervertebral disc

In flexion loading as shown in Figure 5(a), von Mises stress was studied for all the intervertebral disc from C3 to C6 levels, which exhibited a lower range when compared with the simulation results from the literature.³⁴ The C3–C4 level disc exhibited a lower stress compared to the simulation values. The von Mises stress in the disc summary for flexion is given in Figure 3.

In extension loading as shown in Figure 5(b), stress observed on all the segments from C3 to C6 levels were similar to the finite element results of Zafarparandeh et al.³⁴ In the C3–C4 level, the stress across the disc exhibited good agreement with the simulation values, in C4–C5 and C5–C6 levels similar trend was seen.

Figure 5.

Comparison of intervertebral disc stresses for flexion of C3–C6 a) Flexion (b) Extension (c) Right Lateral bending (d) Right Axial rotation.

Disc stress in right lateral bending is shown in Figure 5(c). In the C3–C4 level, the stress across the disc exhibited a higher stress value compared to the published simulation value.³⁴ In the C4–C5 and C5–C6 level discs, the FE analysis exhibited a good agreement and a similar trend was seen with the simulation results.

As shown in Figure 5(d) for right axial rotation, disc stress calculated in C3–C4 level was 2.04 MPa, was higher than the simulation results.³⁴ The maximum disc stress values across C4–C5 level was calculated as 1.9 MPa, was quite higher than the simulation values. In the C5–C6 level, the disc stress exhibited good agreement with the simulation.³⁴ The differences in results with the previously published articles have caused mainly due to the different combination of mechanical properties used in the articles.

Conclusion

In flexion and extension loading from FE studies, it was observed that the C3–C4 level and C5–C6 level exhibit more range of motion in comparison with C4–C5 level. In the lateral bending, all the level from C3 to C6 contribute the same physiological motion by FE studies. In the axial rotation, the range of motion from C3 to C6 levels were observed and it was found that the C4–C5 level is little higher in comparison with the C3–C4 level. In summary, this study presents a validated subject-specific FE Model which ma help clinician in decision making, surgery planning and better treatment for patients.

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

ORCID iD

Sandipan Roy

References

Kurtz

Edidin

Spine technology handbook. New York, NY: Elsevier, 2006.

Hickey

Hukins

DW.

X-ray diffraction studies of the arrangement of collagenous fibres in human fetal intervertebral disc. J Anat 1980; 131(Pt 1): 81.

Yanagishita

Function of proteoglycans in the extracellular matrix. Pathol Int 1993; 43(6): 283–293.

Panjabi

White

AA.

Clinical biomechanics of the spine. Philadelphia, PA: JB Lippincott, 1990

Waheed

MAA

Hasan

Tan

, et al. Cervical spine pathology and treatment: a global overview. J Spine Surg 2020; 6(1): 340.

Raciborski

Gasik

Kłak

Disorders of the spine. A major health and social problem. Reumatologia 2016; 54(4): 196.

Safiri

Kolahi

Hoy

, et al. Global, regional, and national burden of neck pain in the general population, 1990-2017: systematic analysis of the Global Burden of Disease Study 2017. BMJ 2020; 368: m791.

Oxland

TR.

A history of spine biomechanics. Unfallchirurg 2015; 118(Suppl. 1): 80–92.

Erdemir

Guess

Halloran

, et al. Considerations for reporting finite element analysis studies in biomechanics. J Biomech 2012; 45(4): 625–633.

10.

Nikkhoo

Cheng

Wang

, et al. Development and validation of a geometrically personalized finite element model of the lower ligamentous cervical spine for clinical applications. Comput Biol Med 2019; 109: 22–32.

11.

Panjabi

MM.

Cervical spine models for biomechanical research. Spine 1998; 23(24): 2684–2699.

12.

Clavenna

Beutler

Gudipally

, et al. The biomechanical stability of a novel spacer with integrated plate in contiguous two-level and three-level ACDF models: an in vitro cadaveric study. Spine J 2012; 12(2): 157–163.

13.

Li-Ping

Yi-Kai

Bo-Jin

, et al. Morphological changes of the in vitro cervical vertebral canal and its cast form during flexion, extension, and lateral bending. J Manipulative Physiol Ther 2010; 33(2): 132–137.

14.

Kumaresan

Yoganandan

Pintar

, et al. Finite element modeling of the cervical spine: role of intervertebral disc under axial and eccentric loads. Med Eng Phys 1999; 21(10): 689–700.

15.

Anderst

Donaldson

III Lee

, et al. Continuous cervical spine kinematics during in vivo dynamic flexion-extension. Spine J 2014; 14(7): 1221–1227.

16.

Dehner

Schick

Hell

, et al. In-vivo kinematics of the cervical spine in frontal sled tests. Glob J Health Sci 2013; 5(3): 115.

17.

Mao

, et al. Ranges of cervical intervertebral disc deformation during an in vivo dynamic flexion–extension of the neck. J Biomech Eng 2017; 139(6): 0645011–0645017.

18.

Yoganandan

Kumaresan

Voo

, et al. Finite element applications in human cervical spine modeling. Spine 1996; 21(15): 1824–1834.

19.

Panzer

Cronin

DS.

C4–C5 segment finite element model development, validation, and load-sharing investigation. J Biomech 2009; 42(4): 480–490.

20.

Zhang

Teo

HW.

Development and validation of a C0–C7 FE complex for biomechanical study. J Biomech Eng 2005; 127: 729–735.

21.

Kallemeyn

Gandhi

Kode

, et al. Validation of a C2–C7 cervical spine finite element model using specimen-specific flexibility data. Med Eng Phys 2010; 32(5): 482–489.

22.

Zhang

Teo

, et al. Finite element analysis of moment-rotation relationships for human cervical spine. J Biomech 2006; 39(1): 189–193.

23.

Panjabi

Crisco

Vasavada

, et al. Mechanical properties of the human cervical spine as shown by three-dimensional load–displacement curves. Spine 2001; 26(24): 2692–2700.

24.

Musapoor

Nikkhoo

Haghpanahi

A finite element study on intra-operative corrective forces and evaluation of screw density in scoliosis surgeries. Proc Inst Mech Eng Part H: J Eng Med 2018; 232(12): 1245–1254.

25.

Kallemeyn

Tadepalli

Shivanna

, et al. An interactive multiblock approach to meshing the spine. Comput Methods Programs Biomed 2009; 95(3): 227–235.

26.

Biswas

Roy

Pradhan

, et al. Effects of cervical disc replacement and anterior fusion for different bone conditions: a finite element study. Int J Multiscale Comput Eng 2019; 17(4): 411–427.

27.

Bhattacharya

Roy

Rana

, et al. Biomechanical performance of a modified design of dynamic cervical implant compared to conventional ball and socket design of an artificial intervertebral disc implant: a finite element study. J Mech Med Biol 2019; 19(4): 1950017.

28.

Panjabi

Chen

Shin

, et al. The cortical shell architecture of human cervical vertebral bodies. Spine 2001; 26(22): 2478–2484.

29.

Gilad

Nissan

A study of vertebra and disc geometric relations of the human cervical and lumbar spine. Spine 1986; 11(2): 154–157.

30.

Pooni

Hukins

DWL

Harris

, et al. Comparison of the structure of human intervertebral discs in the cervical, thoracic and lumbar regions of the spine. Surg Radiol Anat 1986; 8(3): 175–182.

31.

Pitzen

Schmitz

Georg

, et al. Variation of endplate thickness in the cervical spine. Eur Spine J 2004; 13(3): 235–240.

32.

Yoganandan

Kumaresan

Pintar

FA.

Geometric and mechanical properties of human cervical spine ligaments. J Biomech Eng 2000; 122(6): 623–629.

33.

Arab

Merdji

Benaissa

, et al. Finite-Element analysis of a lateral femoro-tibial impact on the total knee arthroplasty. Comput Methods Programs Biomed 2020; 192: 105446.

34.

Zafarparandeh

Erbulut

Ozer

AF.

Influence of three-dimensional reconstruction method for building a model of the cervical spine on its biomechanical responses: a finite element analysis study. Adv Mech Eng 2016; 8(3): 1687814016638809.

35.

Cai

Sang

Yuchi

, et al. Using finite element analysis to determine effects of the motion loading method on facet joint forces after cervical disc degeneration. Comput Biol Med 2020; 116: 103519.

36.

Traynelis

Donaher

Roach

, et al. Biomechanical comparison of anterior Caspar plate and three-level posterior fixation techniques in a human cadaveric model. J Neurosurg 1993; 79(1): 96–103.

37.

Wheeldon

Pintar

Knowles

, et al. Experimental flexion/extension data corridors for validation of finite element models of the young, normal cervical spine. J Biomech 2006; 39(2): 375–380.

38.

Yoganandan

Pintar

Stemper

, et al. Level-dependent coronal and axial moment-rotation corridors of degeneration-free cervical spines in lateral flexion. JBJS 2007; 89(5): 1066–1074.

39.

Yoganandan

Stemper

Pintar

, et al. Normative segment-specific axial and coronal angulation corridors of subaxial cervical column in axial rotation. Spine 2008; 33(5): 490–496.

40.

SK.

Finite element modeling of multi-level cervical spinal segments (C3–C6) and biomechanical analysis of an elastomer-type prosthetic disc. Med Eng Phys 2006; 28(6): 534–541.