JAMP  Vol.5 No.4 , April 2017
Biomechanical Study of Vertebral Compression Fracture Using Finite Element Analysis
Abstract: This research aimed to mechanically analyze vertebral stress concentration in one healthy subject and one subject with osteoporotic first lumbar (L1) vertebral compression fracture by using finite element analysis (FEA). We constructed three-dimensional image-based finite element (FE) models (Th12L2) by using computed tomographic (CT) digital imaging and communications in medicine (DICOM) for each patient and then conducted exercise stress simulations on the spine models. The loadings on the 12th thoracic vertebra (Th12) due to compression, flexion, extension, lateral bending, and axial rotation were examined within the virtual space for both spine models. The healthy and vertebral compression fracture models were then compared based on the application of equivalent vertebral stress. The comparison showed that vertebral stress concentration increased with all stresses in the vertebral compression fracture models. In particular, compression and axial rotation caused remarkable increases in stress concentration in the vertebral compression fracture models. These results suggest that secondary vertebral compression fractures are caused not only by bone fragility but possibly also by the increase in vertebral stress concentration around the site of the initial fracture

1. Introduction

With the recent aging of the population, the trend toward increasing numbers of patients with osteoporosis has made fragility fractures a major problem in society. Bone fractures cause loss of motor functions and other ailments that markedly lower patient quality of life and lead to invalidity among many of the elderly, making fractures a grave concern in a society which is ageing. The diagnosis of osteoporosis demands methods of discovering and preventing bone fractures early. However, presently, no valid method of diagnosis has been established to quantitatively evaluate the extent of fracture risk due to the mechanical phenomenon of bone destruction. Moreover, the incidence of secondary compression fractures after vertebral fractures is high. This is thought to be caused by bone fragility by osteoporosis. Finite element analysis (FEA) has been introduced in the field of biomechanics over the last few decades. This methodology has been utilized in many clinical applications and gained popularity especially in the prediction of vertebral strength due to its subtle relationships that exist between structure and functionality under a variety of conditions. Moreover, due to the complexity and difficulty of in-vitro and in-vivo experiments, FEA seems to give more promising results. Furthermore, this computational approach reduces the cost and danger of other testing procedures, allowing one to achieve certain individualization when organ geometry and specific loading condition can be fully customized by means of medical image treatment and bio- mechanics simulation technology. The reliability of FEA is subsequently strengthened by the recent finding which demonstrated its better correlations to vertebral strength then Dual-energy X-ray absorptiometry (DXA) approach [1] . In this study, FE models were constructed that reflect the bone density distribution and spinal shape of one healthy subject and one subject with vertebral compression fracture for the purpose of mechanically analyzing vertebral stress concentration. This research aimed to mechanically analyze vertebral stress concentration in one healthy subject and one subject with osteoporotic first lumbar (L1) vertebral compression fracture by using FEA.

2. Materials and Methods

2.1. Patient-Specific FE Modeling

Bone geometrical features were extracted from computed tomographic (CT) di- gital imaging and communications in medicine (DICOM) by using Mechanical Finder (MF) software (Research Center of Computational Mechanics Co. Ltd. Japan) [2] . Individual complex bone shapes and heterogeneous bone density distributions were considered in this bone modeling procedure. Heterogeneous bone density distributions are related to the Young’s modulus of bone and vary between cancellous bones and around the regions between the cortical and cancellous bones. To reflect this heterogeneity in the finite element analysis (FEA), the MF software program was used to calculate the apparent bone density and determine the Young’s modulus of each element separately [3] [4] [5] . Two sets of spinal models of healthy and osteoporotic subjects were developed. Written informed consent, permission, and cooperation were obtained from a 29-year- old Japanese male healthy subject (weight, 78 kg and height, 176 cm) and an 86-year-old Japanese female patient with osteoporosis and first lumbar (L1) vertebral compression fracture (weight, 47 kg and height, 160 cm; fish-type fracture; Figure 1). To create the FE models, CT scan images of the patients’ vertebrae, from the 12th thoracic vertebra (Th12) to the second lumbar vertebra (L2), were obtained. The three-dimensional (3-D) image-based FE models were then constructed based on the extracted bone edges of the region of interest (ROI) around the outer region of the cortical bone on the CT scan images, to obtain the anatomical structure of the spinal bone. Because of the structural complexity of the vertebrae, we adopted tetrahedral elements instead of cubic elements to represent the smooth surface of the spinal bone [6] . The trabecular and inner portion of the cortical bone were modeled by using 3-mm linear tetrahedral elements. Triangular shell elements with a thickness of 0.4 mm were also adopted on the outer surface of the cortex to represent the thin cortical shell [7] . On average, the healthy spinal and osteoporotic models had 804,467 and 790,408 tetrahedral solid elements, and 105,252 and 103,844 triangular shell elements, respectively.

2.2. Calculation of the Bone Material Properties of the Spine FE Models

The bone density of an element was determined from the average CT value (Hounsfield units [HU]) of 17 points, which was composed of the center and four points distributed on four lines connecting the center point to each apex of the tetrahedral element [8] . The bone density of each FE was computed based on the relationship.

Poisson’s ratio was set to a constant value of 0.4, according to that used by Keyak et al. [5] , Reilly and Burstein [9] , and Van Buskirk and Ashman [10] . The elastic modulus of each finite element was determined based on the relationship between Young’s modulus, (Mpa), and the bone density provided by Keyak et al. [11] (Table 1).

Figure 1. Lateral view of plain CT scan images performed during initial examination.

As Young’s modulus is defined by individual elements one by one as described earlier, the heterogeneity of Young’s modulus in the femoral bone can be directly reflected in the FE models. Meanwhile, the yield stress (Mpa) of the models was calculated from the bone density as proposed by Keyak et al. [5] . The final FE models consisted of the T12-L2 vertebrae, intervertebral disks and facet joints (Figure 2). The material properties of the intervertebral disk and facet joint are listed in Table 2.

(a) (b)

Figure 2. Three-dimensional finite-element model: (a) healthy subject (b) osteoporotic subject.

Table 1. Spine material properties.

Table 2. Material properties of the finite-element models.

2.3. Analysis

The FE models were loaded with a compressive force of 1000 N and four rotational/moment loadings on the superior surface of the T12 intervertebral disk to stimulate the four physiological motions/functions of the spine, which represent the movement of flexion, extension, lateral bending, and axial rotation. The inferior side of the L2 intervertebral disk was rigidly fixed. The loading details are listed and depicted in Table 3 and Figure 3, respectively [12] . The biomechanical effects of the osteoporotic bone model were analyzed and compared with those of the healthy bone model. Drucker-Prager stress distributions on the vertebrae were evaluated [13] . In order to evaluate the stress distribution within and between the vertebral bodies, 30 points (10 points for each vertebra) were selected to extract the average Drucker-Prager stress. This point represented a square plate that could measure the average stress distribution distributed uniformly throughout its square volume. The plate was placed in parallel to the vertebral end plates. The distance between each of the plate was set to 5 mm. Figure 4 is Young’s modulus distributions.

The main procedures for FEA statistical analysis are as follows:

1) Bone edges are extracted from the ROI defined from the CT DICOM data.

Figure 3. Five basic vertebral physiological motions.

Figure 4. Young’s modulus distribution.

Table 3. Loading conditions [12] .

2) Three-dimensional image-based FE model contour data are calculated from the bone edge extraction data.

3) Mesh segmentation is applied to the contour data, and a 3-D image-based FE model is constructed.

4) Material properties are bestowed on each element based on conversion equations for the properties of the selected materials.

5) Loads and boundary conditions are established.

6) Statistical processing for each element is conducted based on established methods (i.e., static response analysis and linear regression).

7) From these statistical results, information on aspects such as fracture load, strain distribution, and stress distribution is derived.

3. Results

The load transfer properties (stress and strain) significantly differed between the healthy vertebrae and the osteoporotic vertebrae in five different vertebral physiological motions (Figure 5). In general, the osteoporotic subject tended to pro- duce higher stress and strain than the healthy subject in all physiological movements. The maximum Drucker-Prager stresses (Figure 6) for the healthy subject were 6.45, 10.05, 1.69, 3.03, and 6.48 MPa for compression, flexion, lateral bending, and axial rotation, respectively. Meanwhile, the Drucker-Prager stresses for the osteoporotic subject were 16.85, 10.61, 2.25, 3.90, and 7.96 MPa for compression, flexion, lateral bending, and axial rotation, respectively. The largest relative difference (Figure 7) was found in compression activity (161%), followed by axial rotation (23%), flexion (6%), lateral bending (29%), and extension (33%) activities. The minimum principle strains (Figure 8) for the healthy subject were −3590, −1560, −334, −843, and −637 μs train for compression, flexion, extension, lateral bending, and axial rotation, respectively. Replicating the same distribution pattern as that of the Drucker-Prager stress distributions, the osteoporotic subject had a relatively higher minimum principle strain than the healthy subject. Topping the list (Figure 7) was axial rotation (801%), followed by compression (1816%), flexion (727%), extension (321%), and lateral bending vertebral motion (632%). The result of the average Drucker-Prager stress distributions are shown in Figure 9. The results showed that the greatest Drucker-Prager stress for both subjects was found during compression. For the osteoporotic subject, this stress was substantially higher under relatively similar level of compressive loading, approximately 5 times higher for the osteoporotic than for the healthy subject. It is also important to note that the least relative stress difference was 50% under similar extensive loading, with the osteoporotic subject exhibiting higher stress than the healthy subject. The high degree of these stresses were then correlated with the high

Figure 5. The cross-sectional view of the Drucker-Prager stress distribution on the vertebral body in five different vertebral physiological motions: (a) compression; (b) extension, (c) flexion; (d) lateral bending; and (e) axial rotation.

Figure 6. Drucker-Prager stress on the vertebral body of the healthy and osteoporotic subjects.

degree of principle strain. That is, the Drucker-Prager stress values were directly proportional to the principle strain values, and most of the time, the strains were concentrated in the middle of the trabecular region for each of the vertebrae.

4. Discussion

The recent developments in computational dynamics technology has made possible the use of FEA to conduct a mechanical bone analysis that reflects the complex structural morphology and material properties of bones. Indeed, through

Figure 7. Relative stress difference (%) between the osteoporotic and healthy subjects.

Figure 8. Minimum Principle strain on the vertebral body of the healthy and osteoporotic subjects.

such applications, various mechanical experiments and FEA on vertebral compression fractures have been performed. The method involves conducting FEA structural analysis on 3-D bone models constructed based on DICOM data procured from quantitative CT bone analysis. Then, the intensity of external forces applied in any direction or at any strength are quantitatively examined [14] [15] [16] . First, the CT based 3-D bone structure and bone density distribution can be assessed. Simultaneously, bone density can be evaluated at each tetrahedral element by collecting CT values for bone mass phantoms and plotting them on a calibration curve that correlates the CT value of each tetrahedral element to bone density. Thus, based on the calibration curve, CT values can be converted to bone density data. The bone density data can be arranged just like in a patient’s bone in each area of the 3-D bone model derived from the DICOM data. Computer processing can then reproduce a 3-D bone image with structure and density distribution identical to that of the subject. Further processing can be done to construct FE models, against which virtual loads and boundary

Figure 9. Average Drucker-Prager stress distributions in the vertebral body of the healthy and osteoporotic subjects, measured from the superior end plate of T12 to the inferior end plate of L2: (a) compression, (b) flexion, (c) extension, (d) lateral bending, and (e) axial rotation.

conditions can be established. Then, items such as strain distribution, stress distribution, yield load, and fracture load can be quantitatively assessed. Recently, some attempts have been made to understand the pathology of osteoporosis and judge therapeutic effects in individual patients based on CT data [17] [18] [19] . Other methods have involved the use of implants for surgical evaluations [20] [21] [22] [23] [24] . The first attempt to use FEA in bone strength analysis was made by Brekelmans et al. [25] , whose analysis of the proximal end of the femur was far more accurate than the standard application of beam theory. Keyak [26] published a study in which the accuracy of FEA using CT of the proximal end of the femur was simultaneously tested along with the load testing of vertebrae samples from fresh cadavers. Another study involving lumbar vertebrae samples from fresh cadavers showed that CT-based FEA was highly accurate and correlated highly with vertical compression test results [27] . However, hitherto reports of FE models of osteoporotic vertebrae have mainly focused on individual vertebrae, and no reports have considered vertebrae shape or spinal alignment. In this study, models of a set of three vertebrae were constructed, and further analyses of vertebrae shape and alignment were conducted with respect to vertebral compression fracture models.

Lindsay et al. [28] reported that novel vertebral fractures occur through natural processes at a rate of 19.2% and that after initial vertebrae fractures, the probability of developing secondary compression fractures is extremely high. The results of this study indicated that all the stresses associated with compression, flexion, extension, lateral bending, and axial rotation observed in vertebral compression fracture models resulted in increases in stress concentration in adjacent vertebrae. Thus, it is thought that the occurrence of secondary vertebral compression fractures is due not only to bone fragility but also to alterations in the stress concentrated on adjacent vertebrae due to changes in alignment after the initial vertebral fracture. Considering that hitherto studies on stress loading have not investigated stresses caused by factors besides compression, flexion, and extension [29] [30] and as we have confirmed particularly an increase in stress concentration due to compression and axial rotation, as well as from the perspective of providing guidance on daily activities after vertebral fractures, we believe new insight has been gained.

These phenomena were proven to be well correlated with the deterioration of bone structural strength and reduced bone mass as characterized by osteoporosis. The limitations of this study are a small sample size and we don’t have experimental validations. We will study far more sample after this. Validation study is very difficult, but we hope to try cadaver study.

5. Conclusion

The osteoporotic vertebral model with L1 vertebral compression fracture has significantly affected the load transfer pattern (stress and strain) distributions within the vertebral body. By utilizing the stress and strain distributions of the healthy subject as a comparison tool, we found that the osteoporotic subject seemed to exhibit extremely higher stresses and strains than the healthy subject under the five basic vertebral physiological motions. Worsening this condition was the accompanying uneven stress distribution within and between the vertebral bodies. Therefore, we strongly suggest that for the osteoporotic subject, the risk of vertebral fracture can occur at any time even with daily living activities.


Written informed consent was obtained from the patient for publication of this paper and any accompanying images.

Conflict of Interest

The authors report no conflict of interest concerning the findings specified in this paper.

Cite this paper: Takano, H. , Yonezawa, I. , Todo, M. , Hazli Mazlan, M. , Sato, T. and Kaneko, K. (2017) Biomechanical Study of Vertebral Compression Fracture Using Finite Element Analysis. Journal of Applied Mathematics and Physics, 5, 953-965. doi: 10.4236/jamp.2017.54084.

[1]   Crawford, R.P., Cann, C.E. and Keaveny, T.M. (2003) Finite Element Models Predict in Vitro Vertebral Body Compressive Strength Better than Quantitative Computed Tomography. Bone, 33, 744-750.

[2]   Tawara, D., Sakamoto, J., Murakami, H., Kawahara, N., Oda, J. and Tomita, K. (2010) Mechanical Evaluation by Patient-Specific Finite Element Analysis Demonstrates Therapeutic Effects for Osteoporotic Vertebrae. Journal of the Mechanical Behavior of Biomedical Materials, 3, 31-40.

[3]   Carter, D.R. and Hayes, W.C. (1997) The Compressive Behavior of Bone as Two-Phase of Porous Structure. The Journal of Bone and Joint Surgery, 59, 954-962.

[4]   Keller, T.S. (1994) Predicting the Compressive Mechanical Behavior of Bone. Journal of Biomechanics, 27, 1159-1168.

[5]   Keyak, J.H., Rossia, S.A., Jones, K.A. and Skinner, H.B. (1998) Prediction of Femoral Fracture Load Using Automated Finite Element Modeling. Journal of Biomechanics, 31, 125-133.

[6]   Ulrich, D., Rietbergen, B.V., Weinans, H. and Ruegsegger, P. (1998) Finite Element Analysis of Trabecular Bone Structure: A Comparison of Image-Based Meshing Techniques. Journal of Biomechanics, 31, 1187-1192.

[7]   Imai, K., Ohnishi, I., Yamamoto, S. and Nakamura, K. (2008) In Vivo Assessment of Lumbar Vertebral Strength in Elderly Women Using Computed Tomography-Based Nonlinear Finite Element Model. Spine, 33, 27-32.

[8]   Bessho, M., Ohnishi, I., Matsuyama, J., Matsumoto, T. and Imai, K. (2007) Prediction of Strength and Strain of Proximal Femur by a CT-Based Finite Element Method. Journal of Biomechanics, 3, 31-40.

[9]   Reilly, D.T. and Burstein, A.H. (1975) The Elastic and Ultimate Properties of Compact Bone Tissue. Journal of Biomechanics, 8, 393-405.

[10]   Van Buskirk, W.C.V. and Ashman, R. (1981) The Elastic Moduli of Bone. American Society of Mechanical Engineers, 45, 131-143.

[11]   Keyak, J.H., Lee, I.Y. and Skinner, H.B. (1994) Correlations between Orthogonal Mechanical Properties and Density of Trabecular Bone: Use of Different Densitometric Measures. Journal of Biomedical Materials Research, 28, 1329-1336.

[12]   Han, K.S., Rohlmanna, A., Zandera, T. and Taylorb, W.R. (2013) Lumbar Spinal Loads Vary with Body Height and Weight. Medical Engineering and Physics, 35, 969-977.

[13]   Mazlan, M.H., Todo, M., Takano, H. and Yonezawa, I. (2014) Finite Element Analysis of Osteoporotic Vertebrae with First Lumber (L1) Vertebral Compression Fracture. International Journal of Applied Physics and Mathematics, 4, 267-274.

[14]   Bauera, J.S., Sidorenkoc, I., Muellerd, D., Bauma, T., Isseverb, A.S., Ecksteinf, F., Rummeny, E.J., Link, T.M. and Raeth, C.W. (2014) Prediction of Bone Strength by CT and MDCT-Based Finite-Element Models: How Much Spatial Resolution Is Needed. European Journal of Radiology, 83, 36-42.

[15]   Brouwers, J.E.M., Lambers, F.M., Rietbergen, B.V., Ito, K. and Huiskes, R. (2009) Comparison of Bone Loss Induced by Ovariectomy and Neurectomy in Rats Analyzed by in Vivo Micro-CT. Journal of Orthopaedic Research, 27, 1521-1527.

[16]   Sairyo, K., Goel, V.K., Masuda, A., Vishnubhotla, S., Faizan, A., Biyani, A., Ebraheim, N., Yonekura, D., Murakami, R. and Terai, T. (2006) Three-Dimensional Finite Element Analysis of the Pediatric Lumbar Spine. Part I: Pathomechanism of Apophyseal Bony Ring Fracture. European Spine Journal, 15, 923-929.

[17]   Green, J.O., Diab, T., Allen, M.R., Vidakovic, B., Burr, D.B. and Guldberg, R.E. (2012) Three Years of Alendronate Treatment Does Not Continue to Decrease Microstructural Stresses and Strains Associated with Trabecular Microdamage Initiation beyond Those at 1 Year. Osteoporosis International, 23, 2313-2320.

[18]   Imai, K. (2011) Vertebral Fracture Risk and Alendronate Effects on Osteoporosis Assessed by a Computed Tomography-Based Nonlinear Finite Element Method. Journal of Bone and Mineral Metabolism, 29, 645-651.

[19]   Tawara, D., Sakamoto, J., Murakami, H., Kawahara, N. and Tomita, K. (2011) Patient Specific Finite Element Analyses Detect Significant Mechanical Therapeutic Effects on Osteoporotic Vertebrae during a Three-Year Treatment. Journal of Biomechanical Science and Engineering, 6, 248-261.

[20]   Alizadeh, M., Kadir, M.R.A., Fadhli, M.M., Fallahiarezoodar, A., Azmi, B., Murali, M.R. and Kamarul, T. (2013) The Use of X-Shaped Cross-Link in Posterior Spinal Constructs Improves Stability in Thoracolumbar Burst Fracture: A Finite Element Analysis. Journal of Orthopaedic Research, 31, 1447-1454.

[21]   Choi, K.C., Ryu, K.S., Lee, S.H., Kim, Y.H., Lee, S.J. and Park, C.K. (2013) Biomechanical Comparison of Anterior Lumbar Interbody Fusion: Stand-Alone Interbody Cage versus Interbody Cage with Pedicle Screw Fixation: A Finite Element Analysis. BMC Musculoskeletal Disorders, 14, 220-229.

[22]   Lin, H.M., Pan, Y.N., Liub, C.L., Huangc, L.Y., Huang, C.H. and Chen, C.S. (2013) Biomechanical Comparison of the K-ROD and Dynesys Dynamic Spinal Fixator Systems: A Finite Element Analysis. Bio-Medical Materials and Engineering, 23, 495-505.

[23]   Wilke, H.J., Wenger, K. and Claes, L. (1998) Testing Criteria for Spinal Implants Recommendations for the Standardization of In Vitro Stability Testing of Spinal Implants. European Spine Journal, 7, 148-154.

[24]   Zhang, L., Yang, G., Wub, L. and Yu, B. (2010) The Biomechanical Effects of Osteoporosis Vertebral Augmentation with Cancellous Bone Granules or Bone Cement on Treated and Adjacent Non-Treated Vertebral Bodies: A Finite Element Evaluation. Clinical Biomechanics, 25, 166-172.

[25]   Brekelmans, W.A., Poort, H.W. and Slooff, T.J. (1972) A New Method to Analyze the Mechanical Behavior of Skeletal Parts. Acta Orthopaedica Scandinavica, 43, 301-317.

[26]   Keyak, J.H. (2001) Improved Prediction of Proximal Femoral Fracture Load Using Nonlinear Finite Element Models. Medical Engineering and Physics, 23, 165-173.

[27]   Silva, M.J., Wang, C., Keaveny, T.M. and Hayes, W.C. (1994) Direct and Computed Tomography Thickness Measurements of the Human, Lumbar Vertebral Shell and End Plate. Bone, 15, 409-414.

[28]   Lindsay, R., Silverman, S.L., Cooper, C., Hanley, D.A., Barton, I., Broy, S.B., Licata, A., Benhamou, L., Geusens, P., Flowers, K., Stracke, H. and Seeman, E. (2001) Risk of New Vertebral Fracture in the Year Following a Fracture. Journal of the American Medical Association, 285, 320-323.

[29]   Garo, A., Arnoux, P.J., Wagnac, E. and Aubin, C.E. (2011) Calibration of the Mechanical Properties in a Finite Element Model of a Lumbar Vertebra under Dynamic Compression up to Failure. Medical and Biological Engineering and Computing, 49, 1371-1379.

[30]   Wagnac, E., Arnoux, P.J., Garo, A. and Aubin, C.E. (2012) Finite Element Analysis of the Influence of Loading Rate on a Model of the Full Lumbar Spine under Dynamic Loading Conditions. Medical and Biological Engineering and Computing, 50, 903-915.