Computational Spine

The human spine forms a highly complex system composed of vertebrae, intervertebral discs, ligaments, and muscle groups. Interrogating the function of the spines and the effect of disease on this function is impractical in patients and faces great difficulties in trying to replicate it in the laboratory. The finite element (FE) modeling, allow us to convert the intricate anatomical spinal anatomy structures, obtained via clinical imaging (Computer tomography (CT), magnetic resonance imaging (MRI) and ultrasound (US)), into detailed mathematical models to provide deep insights into the stress distribution, tissue deformation, and the structural integrity of both healthy and compromised spinal segments. In our group, we use FE modeling and classical mechanics of solids, derived from CT and MRI and checked against cadaver experiments, to understand the effect of metastatic spine disease and aging on the function and mechanics of spinal structures, their role in spinal mechanics and to predict the stiffness, strength, and failure of human vertebrae and spinal segments. Our long-term goal is to develop personalized, data-driven, protocols to identify and predict metastatic vertebrae at risk of fracture and the degradation of intervertebral disc mechanical function, to allow clinicians to treat patients prior to neurological compromise and, maintain their daily function and increase their survival. 

This page focuses on the modeling methods. Related experimental and clinical findings appear on our Metastatic Spine [add link when posted], Spine Biomechanics [add link when posted], Intervertebral Disc [add link], and Musculoskeletal Modeling pages.

Metastatic Spine:

A. Homogenized FE model predicted the strength of metastatic vertebrae.

Bone Metastases can destroy bone (osteolytic) via removal of bone and degradation of the bone tissue, and alternatively deposit dense but abnormal bone (osteoblastic), or both (mixed).  Yet most FE models of metastatic vertebrae use none-cancer vertebrae various degree of idealized geometries of osteolytic Foci. Importantly, although breast, prostate and lung primaries produce, representing 2/3 of the population of patients producing spinal bone metastasis, little is known about the effect of either osteoblastic and importantly Mixed lesions  on the mechanism of fracture in pathologic vertebrae and by extension, prediction of this fracture risk. We test whether FE models can predict the strength of vertebrae with real lesions of every type.

In collaboration with professor Philippe Zysset and Dr Marc Stadelman, Zurich University, Switzerland, we showed that when accounting for bone density spatial heterogeneity, application of  conventional homogenized FE model using material properties of healthy bone, predicted the strength of pooled metastatic vertebrae (R² = 0.78) better than bone mineral content, regardless of lesion type [1]. However, lesion quality affected the prediction with the FE prediction explaining 90% and 76% of the strength in mixed and osteolytic vertebrae but only 64% of the strength in osteoblastic vertebrae. Remarkably, 7 of the vertebrae with osteoblastic (Breast and prostate cancer donors) could not be fully tested as their exceeded the testing system upper load limit of 15,000N, exceeding that of the adult human femur, 7,000 N to 13,000 N, considered the strongest bone structure in the body.  A remarkable testament to the effect of bone metastasis on vertebral bone.  Yet in our patient study, 33% of the vertebrae that failed were had osteoblastic or mixed lesions, posing the question how metastatic lesion affect the fracture of human vertebrae?

Homogenized FE prediction of metastatic vertebrae strength

 

 

Figure 1. Homogenized FE model strongly predict the compressive strength of metastatic vertebrae. Note the divergence of vertebrae from breast and lung primaries at the lower strength values and breast and prostate primaries at higher strength from the model, suggesting either additional factors bone material model and or architectural detail) are not fully captured by the FE model. (Stadelmann et al, [1]) 

 

B. Effect of lesion quality on the mechanism of vertebral failure.

To understand how lesion type shapes the way metastatic vertebrae fail, we worked with Raúl Radovitzky's group at MIT to extend a viscoelastic–viscoplastic bone model with continuum damage [2]. Calibrated for each vertebra, the model explained 99% of the variance in measured strength and stiffness across lesion types [1, 2] and revealed marked differences in how osteolytic and osteoblastic vertebrae fail (Figure 2). Osteolytic vertebrae failed through progressive damage largely confined to the lesion, whereas in osteoblastic vertebrae damage spread through low-modulus bone and bypassed the dense regions, which may expose these vertebrae to sudden, burst-type failure. These findings question the clinical view that osteoblastic vertebrae are at low risk of fracture. In both lesion types, failure occurred once damage reached the vertebral cortex, suggesting that cortical involvement should be part of fracture-risk assessment.

Compuational-FE damage model showing effect of lesion on verebral failure

Figure 2 Simulated damage before, at and after failure in (A) an osteolytic and (B) an osteoblastic vertebra. The medial half of each vertebral body is shown transparent. In both, damage begins at the posterior cortex at yield (green arrows, A.I and B.I). In the osteolytic vertebra, it then spreads through a broad region of low-modulus bone (A.II). In the osteoblastic vertebra, it stays within a narrow low-modulus band and largely bypasses the dense, high-modulus bone (B.II). Failure occurs once internal damage joins damage at the anterior and posterior cortex (green and red arrows). (Soltani 2024, Fig. 6)

Our work is now focused on extending the FE model to simulate the effect of lesions on post-failure progression and integrating the intervertebral disc and, lastly, the full spinal unit in the modeling. Predicting where failure begins could, with further validation, help target interventions such as vertebral augmentation to the region at greatest risk.

C. Effect of osteolytic lesion on the functional spinal unit structural response

To better understand the effect of large osteolytic lesions on the state of stress and strain within a spinal unit (defined as two vertebral segments with the intervening disc, representing the basic biomechanical function of the spine), we paired a CT-based FE model of a complete thoracolumbar segment (Figure 3) with strain-gauge measurements on the same specimen (Figure 4), before and after creating a large lytic defect that breached the vertebral cortex [5].

FE model for an FSU model with uncontained defect

Figure 3. FE model of a lumbar spinal unit showing uncontained lesion in the upper vertebra.  (Alkalay. Orthop Res. 2016) 

 

Functional spine unit instrumented with strain gauges

Figure 4  The arrangement of strain gauge rosettes for measuring compressive and shear strains on the vertebral cortex. (Alkalay. Orthop Res. 2016)

The FE model showed that the defect sharply and unevenly raised compressive, tensile, and shear stresses in the bone and cortex around the lesion (Figure 5); tensile (hoop) stress rose by about 350% at the posterior cortex and 250% in the cancellous bone at the lesion, while tensile and shear stresses fell at the endplate of the affected vertebra (Figure 6). The strain gauges confirmed this shift: cortical strains in the lesioned vertebra moved from compression toward tension, consistent with the cortex bulging and beginning to buckle. In the adjacent vertebra, the model predicted little change. However, the gauges recorded higher tensile strains on both sides and uneven compressive strains, showing that the lesion altered load transfer across the disc. The model captured the overall deformation of the segment well but missed many of these local strain changes.

FSU model-figure 5-change in shear strains (-intact-lesion)

Figure 5. The inclusion of the osteolytic lesion markedly increased shear stresses within the lesion vertebra as well as adjacent vertebra.(Alkalay. Orthop Res. 2016)

FSU FE model effect of lesion on stress magnitude

Figure 6 % change in principal and shear stresses magnitude (intact vs. osteolytic) vertebra. (Alkalay. Orthop Res. 2016)

Clinical significance. A vertebra with a large lytic lesion that breaches the cortex may begin to fail under loads of daily living, and the lesion also changes how load reaches the neighboring vertebra, which may raise the risk of adjacent-level collapse. Fracture-risk assessment should therefore account for cortical involvement and the whole spinal segment, not the affected vertebra alone, and FE predictions of local strain need careful experimental validation before clinical use.  Our current work aims to better understand the implication of lesion type  and volumetric geometry on the change in the structural response and risk of failure in metastatic spinal units to allow development of patient-specific risk prediction models.

CT-based structural analysis of spinal strength

Engineering-based mechanics of solids (MOS) models approach,  a simpler alternative to FE, aims to estimate vertebral strength by treating the vertebra as a structural column. Utilizing CT-derived spatial bone mineral density information, this method computes the vertebra's structural axial and bending rigidities, signifying the vertebra's ability to resist applied force and/or bending moment loading within a prescribed deformation limit. Computing these quantities directly from the CT data, the method offers simplicity and rapid computation requiring low computational resources.

Applying this method in non-cancer spinal functional units with simulated lytic defects, this structural analysis strongly predicted vertebral strength under axial compression (r = 0.84–0.85), as compared to a prediction based on bone density or defect size [8]. However, under the more complex axial compression + torsion, the MOS model did less well in predicting the FSU failure loads (r = 0.66) [8].

A significant limitation of the current approach is that the analysis is inherently 2D, constrained to the axial CT image with no structural interaction or influence from directly adjacent (either superior or inferior) bone structures incorporated. Thus, the analysis cannot fully capture BM-mediated spatial changes in vertebral bone structure. Furthermore, the current approach does not include structural indices for shear and torsion, likely explaining the model's low strength prediction under more complex loading conditions, needed to extend this approach to realistic spinal loading. 

Current work:

Figure 7 

Finite element models of the disc

Current FE models of the disc's internal structure, including the boundary between nucleus and annulus and the angle of the annular collagen fibers, are based on studies of healthy discs. Yet both aging and disease markedly alter the disc tissue structure and composition,  which questions the models' ability to faithfully capture and interrogate such changes, thus limiting their use for studying degeneration. In collaboration with the University of Bern (Prof. Philippe Zysset's group), we replaced these assumptions with measurements from MRI [6]. Proton-density images defined the disc's shape and a smooth, element-by-element transition in stiffness from nucleus to annulus, and diffusion tensor imaging (DTI), whose main diffusion direction follows the local collagen fibers, gave every element its own pair of crossing fiber directions. The image-derived fibers reproduced the criss-cross layering of the annulus, with angles of 25° to 50° consistent with histology.

Modeling work flow for FE-DTI

Figure 8. Workflow for converting DTI imaging to FE model for simulation of disc joint mechanical response: Proton-density weighted images [A] are segmented to extract disc tissue material constants, normalized to the highest signal intensity (AF) and the normalized voxel values, λ, assigned to the FE mesh model element [D]. From the DTI data, a k-means clustering algorithm is used to compute and assign the main fiber directions (a4 and a6) [E] to the FE mesh model elements [F]. The FE-DTI model is then used to simulate the in vitro load-displacement curves [G]

A single calibrated model then reproduced the measured response of a human lumbar disc joint in six loading modes (tension, compression, flexion, extension, lateral bending, and axial rotation). Removing the fibers left the model with little stiffness, showing that the DTI-mapped fiber network carries most of the disc's mechanical resistance. 

Our current work aims to further develop this approach via implantation and optimization of the DTI imaging protocol in a clinical MEI scanner to allow imaging of patient discs invivo.   This approach than opens the way to patient-specific disc models that link MRI signs of degeneration to altered spinal mechanics.

Evaluating the effect of Nucleotomy [7]. Nucleotomy is a minimally invasive procedure that removes a portion of the nucleus pulposus  from an injured or herniated intervertebral disc. In collaboration with  Prof. Philippe Zysset's group, our MRI-based disc models  reproduced the 26% increase in motion caused by removing part of the nucleus to within 5% of the measured change in the disc mechanical response [7], demonstrating that anatomically detailed FE models derived from MR imaging can provide insight regrading the effect of surgical procedures  on the resulting changes in the disc mechanical response.  Such information can be utilized for surgical pre-planning to way the cost-benefits of the planned procedure vs. desired post-operation changes in spine mechanical response and stability.

Publications

  1. Stadelmann MA, Schenk DE, Maquer G, Lenherr C, Buck FM, Bosshardt DD, Hoppe S, Theumann N, Alkalay RN, Zysset PK. Conventional finite element models estimate the strength of metastatic human vertebrae despite alterations of the bone's tissue and structure. Bone. 2020;141:115598. doi:10.1016/j.bone.2020.115598
  2. Soltani Z, Xu M, Radovitzky R, Stadelmann MA, Hackney D, Alkalay RN. CT-based finite element simulating spatial bone damage accumulation predicts metastatic human vertebrae strength and stiffness. Front Bioeng Biotechnol. 2024;12:1424553. doi:10.3389/fbioe.2024.1424553
  3. Strack D, Rehtanz N, Soltani Z, Keko M, Subburaj K, Alkalay RN. Effect of CT-based material grouping on finite element strength and stiffness predictions in vertebrae with metastatic lesions. medRxiv. 2026:2026.08.20.26360953. Preprint; under review at Bone. doi:10.64898/2026.08.20.26360953 [Update status when published]
  4. Follet H, Meysen A, Allard V, Gardegaront M, Levillain A, Roux JP, Bermond F, Poncioni S, Schenk D, Stadelmann M, Heidsieck C, Skalli W, Barrey C, Confavreux C, Pialat JB, Alkalay R, Mitton D. Inter-laboratories comparison of different numerical models to four experimental datasets of human vertebrae with or without metastases. J Mech Behav Biomed Mater. 2026;184:107633. doi:10.1016/j.jmbbm.2026.107633
  5. Alkalay RN, Harrigan TP. Mechanical assessment of the effects of metastatic lytic defect on the structural response of human thoracolumbar spine. J Orthop Res. 2016;34(10):1808–1819. doi:10.1002/jor.23154
  6. Stadelmann MA, Maquer G, Voumard B, Grant A, Hackney DB, Vermathen P, Alkalay RN, Zysset PK. Integrating MRI-based geometry, composition and fiber architecture in a finite element model of the human intervertebral disc. J Mech Behav Biomed Mater. 2018;85:37–42. doi:10.1016/j.jmbbm.2018.05.005
  7. Stadelmann MA, Stocker R, Maquer G, Hoppe S, Vermathen P, Alkalay RN, Zysset PK. Finite element models can reproduce the effect of nucleotomy on the multi-axial compliance of human intervertebral discs. Comput Methods Biomech Biomed Engin. 2020;23(13):934–944. doi:10.1080/10255842.2020.1773808
  8. Alkalay R, Adamson R, Miropolsky A, Hackney D. Female human spines with simulated osteolytic defects: CT-based structural analysis of vertebral body strength. Radiology. 2018;288(2):436–444. doi:10.1148/radiol.2018171139