Mechanical loading effects on isthmic spondylolytic lumbar segment: finite element modelling using a personalised geometry

Biomechanics of the isthmic spondylolysis was investigated by using a nonlinear 3D-finite element model (FEM). A personalised in vivo pediatric geometry of L5-S1 low-grade spondylolisthesis patient was used to develop a L5-pelvis motion segment model that took into consideration vertebrae, disc and ligaments. The stress distribution in the affected motion segment under axial force only, and for a combination of flexion and extension was evaluated. Predicted results showed that, under all loading conditions, stresses were much higher on the pedicle and in the dorsal wall of the pars interarticularis due to the abnormal geometry which is consistent with clinical observations.     

Simulation of the behaviour of the L1 vertebra for different material properties and loading conditions

Three-dimensional finite element models of the thoracolumbar junction (T12–L2) and isolated L1 vertebra were developed to investigate the role of material properties and loading conditions on vertebral stresses and strains to predict fracture risk. The geometry of the vertebrae was obtained from computed tomography images. The isolated vertebra model included an L1 vertebra loaded through polymethylmethacrylate plates located at the top and bottom of the vertebra, and the segment model included T12 to L2 vertebrae and seven ligaments, fibrous intervertebral discs and facet joints. Each model was examined with both homogeneous and spatially varying bone tissue properties. Stresses and strains were compared for uniform compression and flexion. Including material heterogeneity remarkably reduced the stiffness of the isolated L1 vertebra and increased the magnitudes of the minimum principal strains and stresses in the mid-transverse section. The stress and strain distributions further changed when physiological loading was applied to the L1 vertebra. In the segment models, including heterogeneous material properties increased the magnitude of the minimum principal strain by 158% in the centre of the mid-transverse section. Overall, the inclusion of heterogeneity and physiological loading increased the magnitude of the strains up to 346% in flexion and 273% in compression.     

Estimation of in vivo inter-vertebral loading during motion using fluoroscopic and magnetic resonance image informed finite element models

Finite element (FE) models driven by medical image data can be used to estimate subject-specific spinal biomechanics. This study aimed to combine magnetic resonance (MR) imaging and quantitative fluoroscopy (QF) in subject-specific FE models of upright standing, flexion and extension. Supine MR images of the lumbar spine were acquired from healthy participants using a 0.5 T MR scanner. Nine 3D quasi-static linear FE models of L3 to L5 were created with an elastic nucleus and orthotropic annulus. QF data was acquired from the same participants who performed trunk flexion to 60° and trunk extension to 20°. The displacements and rotations of the vertebrae were calculated and applied to the FE model. Stresses were averaged across the nucleus region and transformed to the disc co-ordinate system (S1 = mediolateral, S2 = anteroposterior, S3 = axial). In upright standing S3 was predicted to be −0.7 ± 0.6 MPa (L3L4) and −0.6 ± 0.5 MPa (L4L5). S3 increased to −2.0 ± 1.3 MPa (L3L4) and −1.2 ± 0.6 MPa (L4L5) in full flexion and to −1.1 ± 0.8 MPa (L3L4) and −0.7 ± 0.5 MPa (L4L5) in full extension. S1 and S2 followed similar patterns; shear was small apart from S23. Disc stresses correlated to disc orientation and wedging. The results demonstrate that MR and QF data can be combined in a participant-specific FE model to investigate spinal biomechanics in vivo and that predicted stresses are within ranges reported in the literature.     

Search for critical loading condition of the spine–A meta analysis of a nonlinear viscoelastic finite element model

The relative vulnerability of spinal motion segments to different loading combinations remains unknown. The meta-analysis described here using the results of a validated L2–L3 nonlinear viscoelastic finite element model was designed to investigate the critical loading and its effect on the internal mechanics of the human lumbar spine. A Box-Behnken experimental design was used to design the magnitude of seven independent variables associated with loads, rotations and velocity of motion. Subsequently, an optimization method was used to find the primary and secondary variables that influence spine mechanical output related to facet forces, disc pressure, ligament forces, annulus matrix compressive/shear stresses and anulus fibers strain. The mechanical responses with respect to the two most-relevant variables were then regressed linearly using the response surface quadratic model. Axial force and sagittal rotation were identified as the most-relevant variables for mechanical responses. The procedure developed can be used to find the critical loading for finite element models with multi input variables. The derived meta-models can be used to predict the risk associated with various loading parameters and in setting safer load limits.     

Predictive modelling of cervical disc implant wear

This study presents a chain of simulations aimed at estimating the wear in a cervical disc implant and providing insight into the in vivo biomechanical performance of the implant. The simulation chain can start with determining a representative maximum range of motion (ROM) of a person's head. The ROM is used as motion input to a kinematic simulation of the cervical spine containing a disc implant. The cervical spine geometry is obtained from computed tomography (CT) scans and converted to STL format using reverse engineering software. The time histories of the loads imposed by the adjacent vertebrae on the implant, as well as the vertebral relative rotations can be extracted from the kinematic simulation. Alternatively, force and motion profiles prescribed by wear test protocols (e.g. ISO 18192-1 and ASTM F2423-05) can be used. The force and motion profiles are applied as boundary conditions to a non-linear finite element model (FEM) of the implant to determine the time-varying contact stress and slip velocity distributions at the interface between the two halves of the implant. The stresses and slip velocities are used in a linear wear model to estimate the wear rate distribution at the FEM's nodal points where contact occurs. Reverse engineering software is used to triangulate the contact surface so that the total wear volume can be calculated. The simulation chain's predicted wear rate shows good agreement with in vitro results in the literature. The simulation chain is thereby demonstrated to be suitable for comparative pre-experimental studies of spinal implant designs.     

Subject-specific finite element models of long bones: An in vitro evaluation of the overall accuracy

The determination of the mechanical stresses induced in human bones is of great importance in both research and clinical practice. Since the stresses in bones cannot be measured non-invasively in vivo, the only way to estimate them is through subject-specific finite element modelling. Several methods exist for the automatic generation of these models from CT data, but before bringing them in the clinical practice it is necessary to assess their accuracy in the predictions of the bone stresses. Particular attention should be paid to those regions, like the epiphyseal and metaphyseal parts of long bones, where the automatic methods are typically less accurate. Aim of the present study was to implement a general procedure to automatically generate subject-specific finite element models of bones from CT data and estimate the accuracy of this general procedure by applying it to one real femur. This femur was tested in vitro under five different loading scenarios and the results of these tests were used to verify how the adoption of a simplified two-material homogeneous model would change the accuracy with respect to the density-based inhomogeneous one, with special attention paid to the epiphyseal and metaphyseal proximal regions of the bone. The results showed that the density-based inhomogeneous model predicts with a very good accuracy the measured stresses (R(2)=0.91, RMSE=8.6%, peak error=27%), while the two-material model was less accurate (R(2)=0.89, RMSE=9.6%, peak error=35%). The results showed that it is possible to automatically generate accurate finite element models of bones from CT data and that the strategy of material properties mapping has a significant influence on its accuracy.     

Effect of disc degeneration on the muscle recruitment pattern in upright posture: a computational analysis

Based on the sensor driving control mechanism model, the effect of disc degeneration on the trunk muscle recruitment (TMR) pattern was analysed in erect standing posture. A previously developed computational model was used for this analysis, with modifications incorporating the T12-L1 motion segment and additional muscle fascicles. To generate disc degeneration at three different levels (L3–L4, L4–L5, or L5–S1), the material properties of the ground matrix of the annulus and bulk modulus of the nucleus were reduced. The finite element method combined with an optimization technique was applied to calculate the muscle forces. Minimization of deviations in the averaged tensile stress in the annulus fibres at the outermost layer in the five discs was selected for muscle force calculations. The results indicated that the disc degeneration noticeably increased the activation of the superficial muscle (IT and R) even though there was no clear change in the longissimus thoracis. Unlike some of the superficial muscles, activation in the deep muscles (multifidus (ML, MS, MT), LL and Q) was decreased. The change in TMR pattern generated an intervertebral disc angle difference and nucleus pressure increased in the upper level. These differences are expected to be functional in that they reduce the stress at the degenerated disc by changing the muscle activation, which slows down the progress of disc degeneration.     

In-vivo motion characteristics of lumbar vertebrae in sagittal and transverse planes

Lumbar vertebrae are complicated in structure and function. The purpose of this study was to investigate the in-vivo motion characteristics of different portions of the lumbar vertebrae during functional activities. Motion of L2, L3 and L4 was reproduced using a combined dual fluoroscopic and MR imaging technique during flexion–extension and left–right twisting of the trunk. The ranges of motion (ROM) of the proximal vertebra with respect to the distal one at 3 representative locations: the center of the vertebral body, the center of the spinal canal and the tip of the spinous process were measured. Centers of rotation (COR) of the vertebrae were then determined by calculation of the points of zero motion in 2D sagittal and transverse planes. During flexion–extension, the center of the vertebral body moved less than 0.6 mm, while the tip of the spinous process moved less than 7.5 mm in the sagittal plane. The CORs of both L23 (L2 with respect to L3) and L34 were located inside the vertebral body, at a distance about one-third the length of the vertebral body from the posterior edge. During left–right twisting, the center of the vertebral body moved less than 1.0 mm, while the tip of the spinous process moved less than 1.6 mm in the transverse plane. The CORs of both L23 and L34 were located approximately 30 mm anterior to the front edge of the vertebral body. The results of this study may be used to define the ideal locations for surgical placement of the disc prosthesis, thus help improve the prosthesis design and surgical treatment of various pathological conditions.     

The effect of three-dimensional geometrical changes during adolescent growth on the biomechanics of a spinal motion segment

During adolescent growth, vertebrae and intervertebral discs undergo various geometrical changes. Although such changes in geometry are well known, their effects on spinal stiffness remains poorly understood. However, this understanding is essential in the treatment of spinal abnormalities during growth, such as scoliosis.A finite element model of an L3–L4 motion segment was developed, validated and applied to study the quantitative effects of changing geometry during adolescent growth on spinal stiffness in flexion, extension, lateral bending and axial rotation. Height, width and depth of the vertebrae and intervertebral disc were varied, as were the width of the transverse processes, the length of the spinous process, the size of the nucleus, facet joint areas and ligament size. These variations were based on average growth data for girls, as reported in literature.Overall, adolescent growth increases the stiffness with 36% (lateral bending and extension) to 44% (flexion). Two thirds of this increase occurs between 10 and 14 years of age and the last third between 14 years of age and maturity.Although the height is the largest geometrical change during adolescent growth, its effect on the biomechanics is small. The depth increase of the disc and vertebrae significantly affects the stiffness in all directions, while the width increase mainly affects the lateral bending stiffness. Hence, when analysing the biomechanics of the growing adolescent spine (for instance in scoliosis research), the inclusion of depth and width changes, in addition to the usually implemented height change, is essential.     

Finite element lumbar spine facet contact parameter predictions are affected by the cartilage thickness distribution and initial joint gap size

Current finite element modeling techniques utilize geometrically inaccurate cartilage distribution representations in the lumbar spine. We hypothesize that this shortcoming severely limits the predictive fidelity of these simulations. Specifically, it is unclear how these anatomically inaccurate cartilage representations alter range of motion and facet contact predictions. In the current study, cadaveric vertebrae were serially sectioned, and images were taken of each slice in order to identify the osteochondral interface and the articulating surface. A series of custom-written algorithms were utilized in order to quantify each facet joint's three-dimensional cartilage distribution using a previously developed methodology. These vertebrae-dependent thickness cartilage distributions were implemented on an L1 through L5 lumbar spine finite element model. Moments were applied in three principal planes of motion, and range of motion and facet contact predictions from the variable thickness and constant thickness distribution models were determined. Initial facet gap thickness dimensions were also parameterized. The data indicate that the mean and maximum cartilage thickness increased inferiorly from L1 to L5, with an overall mean thickness value of 0.57 mm. Cartilage distribution and initial facet joint gap thickness had little influence on the lumbar range of motion in any direction, whereas the mean contact pressure, total contact force, and total contact area predictions were altered considerably. The data indicate that range of motion predictions alone are insufficient to establish model validation intended to predict mechanical contact parameters. These data also emphasize the need for the careful consideration of the initial facet joint gap thickness with respect to the spinal condition being studied.     

Search for critical loading condition of the spine--a meta analysis of a nonlinear viscoelastic finite element model

The relative vulnerability of spinal motion segments to different loading combinations remains unknown. The meta-analysis described here using the results of a validated L2-L3 nonlinear viscoelastic finite element model was designed to investigate the critical loading and its effect on the internal mechanics of the human lumbar spine. A Box-Behnken experimental design was used to design the magnitude of seven independent variables associated with loads, rotations and velocity of motion. Subsequently, an optimization method was used to find the primary and secondary variables that influence spine mechanical output related to facet forces, disc pressure, ligament forces, annulus matrix compressive/shear stresses and anulus fibers strain. The mechanical responses with respect to the two most-relevant variables were then regressed linearly using the response surface quadratic model. Axial force and sagittal rotation were identified as the most-relevant variables for mechanical responses. The procedure developed can be used to find the critical loading for finite element models with multi input variables. The derived meta-models can be used to predict the risk associated with various loading parameters and in setting safer load limits.     

Finite element estimates of interface stress in the trans-tibial prosthesis using gap elements are different from those using automated contact

When compared with automated contact methods of finite element (FE) analyses, gap elements have certain inherent disadvantages in simulating large slip of compliant materials on stiff surfaces. However, automated contact has found limited use in the biomechanical literature. A non-linear, three-dimensional, geometrically accurate, FE analysis of the trans-tibial limb-socket prosthetic system was used to compare an automated contact interface model with a gap element model, and to evaluate the sensitivity of automated contact to interfacial coefficient of friction (COF). Peak normal stresses and resultant shear stresses were higher in the gap element model than in the automated contact model, while the maximum axial slip was less. Under proximally directed load, compared with automated contact, gap elements predicted larger areas of stress concentration that were located more distally. Gap elements did not predict any relative slip at the distal end, and also transmitted a larger proportion of axial load as shear stress. Both models demonstrated non -linear sensitivity to COF, with larger variation at lower magnitudes of COF. By imposing physical connections between interface surfaces, gap elements distort the interface stress distributions under large slip. Automated contact methods offer an attractive alternative in applications such as prosthetic FE modeling, where the initial position of the limb in the socket is not known, where local geometric features have high design significance, and where large slip occurs under load.     

Biomechanical assessment of stability in the metastatic spine following percutaneous vertebroplasty: effects of cement distribution patterns and volume

Percutaneous vertebroplasty is a minimally invasive, radiologically guided procedure whereby bone cement is injected into structurally weakened vertebrae to provide added biomechanical stability. In addition to treating osteoporotic vertebral fractures, this technique is also used to relieve pain by stabilizing metastatically compromised vertebrae that are at risk of pathologic burst fracture. Optimal cement distribution patterns to improve biomechanical stability to metastatically involved vertebral bodies remain unknown. This study aimed to determine the effect of cement location and volume of cement injected during percutaneous vertebroplasty on improving vertebral stability in a metastatically-compromised spinal motion segment using a parametric poroelastic finite element model. A three-dimensional parametric finite element model of a thoracic spinal motion segment was developed and analyzed using commercially available software. A total of 16 metastatic pre and post vertebroplasty scenarios were investigated using a serrated spherical representation of tumor tissue and various geometric representations of polymethylmethacrylate (PMMA). The effect of vertebroplasty on vertebral bulge, a measure of posterior vertebral body wall motion as an indicator of burst fracture initiation, was assessed. In all cases, vertebroplasty reduced vertebral bulge, but the risk of the initiation of burst fracture was minimized with cement located posterior to the tumor, near the posterior vertebral body wall. Vertebral bulge decreased by up to 62% with 20% cement injection. These findings demonstrate that location and distribution of cement within the vertebral body has a noticeable effect on the restoration of biomechanical stability following percutaneous vertebroplasty.     

Finite Element Analysis of a New Pedicle Screw-Plate System for Minimally Invasive Transforaminal Lumbar Interbody Fusion

Purpose

Minimally invasive transforaminal lumbar interbody fusion (MI-TLIF) is increasingly popular for the surgical treatment of degenerative lumbar disc diseases. The constructs intended for segmental stability are varied in MI-TLIF. We adopted finite element (FE) analysis to compare the stability after different construct fixations using interbody cage with posterior pedicle screw-rod or pedicle screw-plate instrumentation system.

Methods

A L3–S1 FE model was modified to simulate decompression and fusion at L4–L5 segment. Fixation modes included unilateral plate (UP), unilateral rod (UR), bilateral plate (BP), bilateral rod (BR) and UP+UR fixation. The inferior surface of the S1 vertebra remained immobilized throughout the load simulation, and a bending moment of 7.5 Nm with 400N pre-load was applied on the L3 vertebra to recreate flexion, extension, lateral bending, and axial rotation. Range of motion (ROM) and Von Mises stress were evaluated for intact and instrumentation models in all loading planes.

Results

All reconstructive conditions displayed decreased motion at L4–L5. The pedicle screw-plate system offered equal ROM to pedicle screw-rod system in unilateral or bilateral fixation modes respectively. Pedicle screw stresses for plate system were 2.2 times greater than those for rod system in left lateral bending under unilateral fixation. Stresses for plate were 3.1 times greater than those for rod in right axial rotation under bilateral fixation. Stresses on intervertebral graft for plate system were similar to rod system in unilateral and bilateral fixation modes respectively. Increased ROM and posterior instrumentation stresses were observed in all loading modes with unilateral fixation compared with bilateral fixation in both systems.

Conclusions

Transforaminal lumbar interbody fusion augmentation with pedicle screw-plate system fixation increases fusion construct stability equally to the pedicle screw-rod system. Increased posterior instrumentation stresses are observed in all loading modes with plate fixation, and bilateral fixation could reduce stress concentration.     

Biomechanical analysis of lumbar interbody fusion cages with various lordotic angles: a finite element study

Inappropriate lordotic angle of lumbar fusion cage could be associated with cage damage or subsidence. The biomechanical influence of cage lordotic angle on lumbar spine has not been fully investigated. Four surgical finite element models were constructed by inserting cages with various lordotic angles at L3-L4 disc space. The four motion modes were simulated. The range of motion (ROM) decreased with increased lordotic angle of cage in flexion, extension, and rotation, whereas it was not substantially changed in bending. The maximum stress in cage decreased with increased lordotic angle of cage in all motion modes. The maximum stress in endplate at surgical level increased with increased lordotic angle of cage in flexion and rotation, whereas it was not substantially changed in extension and bending. The facet joint force (FJF) was much smaller than that for the intact conditions in extension, bending, and rotation, while it was not substantially changed in flexion. In conclusion, the ROM, stresses in the cage and endplate at surgical level are sensitive to the lordotic angle of cage. The increased cage lordotic angle may provide better stability and reduce the risk of cage damage, whereas it may increase the risk of subsidence in flexion and rotation.     

Possible role of stresses in inducing spinal stenosis--a long term complication following disk excision

A three-dimensional finite element model of an intact ligamentous lumbar motion segment (L3-4) was used to predict stresses in the pars interarticularis regions of the modeled vertebral bodies. The changes in stresses following disk excision, as compared to the intact model, also were computed. The predicted results show an increase in stresses in the posterior bony elements following disk excision. In some patients over a long period of time this increase in stresses, in association with other clinical factors, may lead to bony hypertrophy of the structures that surround the nerve roots. Ultimately, over a long period of time the increase in pressure on the entrapped nerve root may induce recurrent pain and other complications reported in the literature.     

Effect of intervertebral translational flexibilities on estimations of trunk muscle forces,kinematics, loads,and stability

Due to the complexity of the human spinal motion segments, the intervertebral joints are often simulated in the musculoskeletal trunk models as pivots thus allowing no translational degrees of freedom (DOFs). This work aims to investigate, for the first time, the effect of such widely used assumption on trunk muscle forces, spinal loads, kinematics, and stability during a number of static activities. To address this, the shear deformable beam elements used in our nonlinear finite element (OFE) musculoskeletal model of the trunk were either substantially stiffened in translational directions (SFE model) or replaced by hinge joints interconnected through rotational springs (HFE model). Results indicated that ignoring intervertebral translational DOFs had in general low to moderate impact on model predictions. Compared with the OFE model, the SFE and HFE models predicted generally larger L4–L5 and L5–S1 compression and shear loads, especially for tasks with greater trunk angles; differences reached ～15% for the L4–L5 compression, ～36% for the L4–L5 shear and ～18% for the L5–S1 shear loads. Such differences increased, as location of the hinge joints in the HFE model moved from the mid-disc height to either the lower or upper endplates. Stability analyses of these models for some select activities revealed small changes in predicted margin of stability. Model studies dealing exclusively with the estimation of spinal loads and/or stability may, hence with small loss of accuracy, neglect intervertebral translational DOFs at smaller trunk flexion angles for the sake of computational simplicity.     

Compression of interventricular septum during right ventricular pressure loading

The interventricular septum, which flattens and inverts in conditions such as pulmonary hypertension, is considered by many to be an unstressed membrane, in that its position is assumed to be determined solely by the transseptal pressure gradient. A two-dimensional finite element model was developed to investigate whether compression and bending moments (behavior incompatible with a membrane) exist in the septum during diastole under abnormal loading, i.e., pulmonary artery (PA) constriction. Hemodynamic and echocardiographic data were obtained in six open-chest anesthetized dogs. For both control and PA constriction, the measured left ventricular and right ventricular pressures were applied to a residually stressed mesh. Adjustments were made to the stiffness and end-bending moments until the deformed and loaded residually stressed mesh matched the observed configuration of the septum. During PA constriction, end-bending moments were required to obtain satisfactory matches but not during control. Furthermore, substantial circumferential compressive stresses developed during PA constriction. Such stresses might impede septal blood flow and provoke the unexplained ischemia observed in some conditions characterized by abnormal septal motion.     

The effects of different articulate curvature of artificial disc on loading distribution

Purpose: Deeper insights into the mechanical behavior of lumbar disc prostheses are required. Prior studies on the biomechanical performance of artificial discs were mostly performed with finite element analyses, but this has never been analyzed with altering articulate curvature. This study aimed to ascertain the influence of the geometry of a ball-and-socket disc prosthesis for the lumbar spine. Materials and Methods: Three-dimensional finite element model of human L4-L5 was reconstructed. Convex, concave, and elliptic artificial disc models were also established with Computer-Aided-Design software. Simulations included: (1) three articulate types of polyethylene (PE) insert were implanted inferiorly and (2) concave and convex PE inserts were implanted on the superior or inferior sides in flexion/extension, lateral bending, and axial rotation in the lumbar spine. Shear stresses and von Mises stresses on PE insert were assessed for their loading distributions. Results: High shear stresses of all articulate types occurred in flexion, and convex PE insert performed the maximum stress of 23.81 MPa. Under all conditions, stresses on concave PE inserts were distributed more evenly and lower than those on the convex type. Elliptic geometry enabled confining the rotation of the motion unit. The shear force on the convex PE insert on the inferior side could induce transverse crack because the shear stress exceeded yielding shear stress. Conclusions: The concave PE insert on the inferior side not only decreased loading concentration but had relatively low stress. Such a design may be applicable for artificial discs.     

Effect of a Double Helical Spring Decompression Structure Backpack on the Lumbar Spine Biomechanics of School-Age Children: A Finite Element Study

Background: A children’s backpack is one of the important school supplies for school-age children. Long-term excessive weight can cause spinal deformity that cannot be reversed. This study compared a double helical spring decompression structure backpack (DHSB) with a traditional backpack (TB) to explore the optimization of decompression devices on upper body pressure. The finite element (FE) method was then used to explore the simulation of lumbar stress with different backpacks, in order to prove that DHSB can reduce the influence of backpack weight on lumbar vertebrae, avoid the occurrence of muscle discomfort and spinal deformity in children; Methods: 18 male children subjects (age: 12.5 ± 0.6 years; height: 145.5 ± 1.9 cm; bodyweight: 40.8 ± 3.1 kg) ran with DHSB and TB at a speed of 3.3 ± 0.2 m/s. Flexible pressure sensors were used to measure the pressure on the shoulder, back, and waist during running. The pressure data was then inputted into the FE model to simulate the effect of carrying different backpacks on the stress of the lumbar intervertebral disc (IVD); Result: There was a significant difference in shoulder and waist peak pressure between the DHSB and TB during the running posture. At a speed of 3.3 ± 0.2 m/s, the peak pressure of the shoulder and waist decreased. After finite element analysis, it was found that carrying DHSB on the back could effectively reduce the intervertebral disc pressure between L4-5 and L5-S1 by 27.9% and 34.1%, respectively; Conclusion: DHSB can effectively reduce the pressure on the shoulder and waist when children are running and can reduce the influence of backpacks on children’s posture to a certain extent. By finite element analysis, it is found that carrying DHSB can effectively reduce the stress of the lumbar intervertebral disc, and the damage to lumbar vertebrae is lower than with a TB.