Biomechanical Analysis of the Spine in Diffuse Idiopathic Skeletal Hyperostosis: Finite Element Analysis

: Patients with diffuse idiopathic skeletal hyperostosis (DISH) develop fractures of the vertebral bodies, even in minor trauma, because of the loss of ﬂexibility, which causes difﬁculties in fusing vertebrae; therefore, the diagnosis of spine injuries may be delayed. We used the three-dimensional ﬁnite element method to add data on ossiﬁcation to the healthy vertebral model in order to investigate how stress in intervertebral discs changes with bone shape and whether these changes present any risk factors. A healthy spine model and a DISH ﬂat model (T8–sacrum) were generated from medical images. As an ossiﬁed hypertrophic model, T11–T12 was cross-linked with hypertrophic ossiﬁcation, and hypertrophy was found to be 5 and 10 mm. An ossifying hypertrophic groove model (5 and 10 mm) was created at T11–T12 and T11–L1. A groove was created at the center of T12, and the radius of curvature of the groove was set to 1 and 2.5 mm. An extension force and ﬂexion force were applied to the upper part of T8, assuming that external forces in the direction of ﬂexion and extension were applied to the spine. Stresses were greater in the DISH ﬂat model than in the healthy model. In the hypertrophic ossiﬁcation model, the stress on the vertebral body was similar to greater ossiﬁcation in extension and ﬂexion. In the ossiﬁed hypertrophic groove model, the stress at the center of the groove increased. In DISH, vertebrae are more susceptible to stress. Furthermore, depending on the morphology of ossiﬁcation, stresses on the vertebrae and intervertebral discs differed even with similar loads. An examination of ossiﬁcation geometry may help surgeons decide the thoracolumbar spine’s stress elevated position in patients with DISH, thereby contributing to the understanding of the pathogenesis of pain.


Introduction
Diffuse idiopathic skeletal hyperostosis (DISH) was defined by Resnick et al. in 1975 [1]. The most commonly used diagnostic criteria are the involvement of at least four contiguous vertebrae of the thoracic spine, the preservation of intervertebral disc spaces, and the absence of gross degeneration or fusion of the apophyseal and sacroiliac joints. DISH most commonly affects the spine and often presents as back pain and stiffness. Furthermore, it is associated with fractures of the vertebral body, even with minor trauma, because of the loss of flexibility, which is problematic [2][3][4]. Computed tomography (CT) has made it possible to easily diagnose the morphology of ossification [5,6]. Previous studies added the thickness of the bony ridge as a diagnostic criterion [7]. However, the morphology of ossification in practice may be a combination of ventral overhang and flatness ( Figure 1). There is currently no classification for the various forms of ossification; however, we hypothesized that under similar loads and traumatic maneuvers, the morphology of ossification and the depth of the ossification groove may vary and also that the stresses applied to the vertebrae may differ. A stress analysis using the ossification geometry of DISH has not yet been conducted. Therefore, we used the three-dimensional finite element method (3D-FEM) to add ossification to the same vertebral model in order to investigate how the stress on intervertebral discs changes with bone shape and whether these changes present any risk factors.
ppl. Sci. 2021, 11, x FOR PEER REVIEW 2 of 13 ( Figure 1). There is currently no classification for the various forms of ossification; however, we hypothesized that under similar loads and traumatic maneuvers, the morphology of ossification and the depth of the ossification groove may vary and also that the stresses applied to the vertebrae may differ. A stress analysis using the ossification geometry of DISH has not yet been conducted. Therefore, we used the three-dimensional finite element method (3D-FEM) to add ossification to the same vertebral model in order to investigate how the stress on intervertebral discs changes with bone shape and whether these changes present any risk factors.

Patient Images
CT images (slice thickness of 0.6 mm) of the spine, from the thoracic spine to the pelvis, of an adult man were obtained with the Brilliance 64 CT scanner (Philips Healthcare, Amsterdam, the Netherlands). The present study was approved by the Ethics Committee at the Center for Clinical Research, Yamaguchi University Hospital (Ube, Japan; approval No. H28-054). Written informed consent for this study and its publication was obtained from all subjects.

FEM Model Construction
Model construction was performed using finite element analysis software (Simpleware ScanIP, version M-2017.06; Synopsys Inc., Mountain View, CA, USA). After the spine was extracted, vertebrae were mapped as cancellous and cortical bones, and the sizes of the intervertebral discs were adjusted to match the sizes of the endplates of each vertebra. When cortical and cancellous bones were being separated, cortical bone was delineated using the CT values (>1000 Hounsfield units) of cortical bone. Cancellous bone was colorcoded on the assumption that it was contained within cortical bone. Due to the small size of the facet joints, the computer was unable to automatically separate them. Therefore, manual distinctions were performed when checking CT scans. A 3D whole-spine model was constructed by individually mapping all vertebrae and intervertebral discs from the thoracic to sacral regions. The gap between each vertebra and intervertebral disc was considered to be completely restricted in movement. Facet-joint spaces were created at all levels to allow each vertebra to move independently. The anterior longitudinal ligament,

Patient Images
CT images (slice thickness of 0.6 mm) of the spine, from the thoracic spine to the pelvis, of an adult man were obtained with the Brilliance 64 CT scanner (Philips Healthcare, Amsterdam, the Netherlands). The present study was approved by the Ethics Committee at the Center for Clinical Research, Yamaguchi University Hospital (Ube, Japan; approval No. H28-054). Written informed consent for this study and its publication was obtained from all subjects.

FEM Model Construction
Model construction was performed using finite element analysis software (Simpleware ScanIP, version M-2017.06; Synopsys Inc., Mountain View, CA, USA). After the spine was extracted, vertebrae were mapped as cancellous and cortical bones, and the sizes of the intervertebral discs were adjusted to match the sizes of the endplates of each vertebra. When cortical and cancellous bones were being separated, cortical bone was delineated using the CT values (>1000 Hounsfield units) of cortical bone. Cancellous bone was colorcoded on the assumption that it was contained within cortical bone. Due to the small size of the facet joints, the computer was unable to automatically separate them. Therefore, manual distinctions were performed when checking CT scans. A 3D whole-spine model was constructed by individually mapping all vertebrae and intervertebral discs from the thoracic to sacral regions. The gap between each vertebra and intervertebral disc was considered to be completely restricted in movement. Facet-joint spaces were created at all levels to allow each vertebra to move independently. The anterior longitudinal ligament, posterior longitudinal ligament, and ligamentum flavum were also added to the model ( Figure 2). The contact of the intervertebral disc, vertebral body, ligament, and vertebral body/intervertebral disc were used to make translational and rotational motion as well as all other active degrees of freedom equal for a pair of surfaces. Moreover, all facet joints were independently constructed. Each facet joint was separated from the other facet, and the surface was set to have a coefficient of friction of 0 in order to account for joint fluid in case of contact. This model was used as the healthy spine model. In this analysis, the cranial end was T8 and the caudal side was the sacrum. All elements were considered to be linear elastic materials. Young's modulus was used as a reference, as previously described [8][9][10][11][12][13] (Table 1).  [9], Cowin [10], Pitzen et al. [9], and Ottardi et al. [13]. The contact of the intervertebral disc, vertebral body, ligament, and vertebral body/ intervertebral disc were used to make translational and rotational motion as well as all other active degrees of freedom equal for a pair of surfaces. Moreover, all facet joints were independently constructed. Each facet joint was separated from the other facet, and the surface was set to have a coefficient of friction of 0 in order to account for joint fluid in case of contact. This model was used as the healthy spine model. In this analysis, the cranial end was T8 and the caudal side was the sacrum. All elements were considered to be linear elastic materials. Young's modulus was used as a reference, as previously described [8][9][10][11][12][13] ( Table 1).

Model for DISH
In the DISH flat model, the physical properties of the whole anterior longitudinal ligament were changed to Young's modulus of cortical bone ( Figure 3A).

Model for DISH
In the DISH flat model, the physical properties of the whole anterior longitudinal ligament were changed to Young's modulus of cortical bone ( Figure 3A). As an ossified hypertrophic model, T11-T12 was cross-linked with hypertrophic ossification, and hypertrophy was modeled with a thickness of 5 mm (T-5) ( Figure 3b) and 10 mm (T-10) (Figure 3c). In addition, an ossifying hypertrophic groove model (5 and 10 mm) was created at T11-T12 and T12-L1. A groove was created at the center of T12, and the radius of curvature of this groove was set to 1 and 2.5 mm (D-5/R-2.5, D-10/R-2.5, D-5/R-1, D-10/R-1, respectively) ( Figure 4a,b).
In the healthy spine model, the total numbers of elements and nodes were 903,621 and 203,312 using tetrahedral elements, respectively. The mesh in the groove area was small to allow the analysis to converge ( Table 2).   As an ossified hypertrophic model, T11-T12 was cross-linked with hypertrophic ossification, and hypertrophy was modeled with a thickness of 5 mm (T-5) ( Figure 3b) and 10 mm (T-10) (Figure 3c). In addition, an ossifying hypertrophic groove model (5 and 10 mm) was created at T11-T12 and T12-L1. A groove was created at the center of T12, and the radius of curvature of this groove was set to 1 and 2.  . Ossifying hypertrophic grooves at 5 (a) and 10 mm (b) were created at T11-T12 and T12-L1. In these diffuse idiopathic skeletal hyperostosis models, the sacrum was fixed and the cranial endplate of T8 was displaced in the direction of flexion and extension.

Load Application
With the sacrum completely restrained in all directions, a forced displacement of 30 mm to ventral and dorsal direction was applied to the cranial endplate of T8, assuming extension and flexion ( Figure 4b).
An analysis was performed using Patran and MARC (MSC Software, Newport Beach, CA, USA). Sixteen different compression combinations were evaluated, a mechanical analysis of vertebral bodies and discs was performed based on the visual findings, and the highest value for maximum von Mises stress was recorded for each combination. The unit is the megapascal (MPa). Ossifying hypertrophic grooves at 5 (a) and 10 mm (b) were created at T11-T12 and T12-L1. In these diffuse idiopathic skeletal hyperostosis models, the sacrum was fixed and the cranial endplate of T8 was displaced in the direction of flexion and extension.
In the healthy spine model, the total numbers of elements and nodes were 903,621 and 203,312 using tetrahedral elements, respectively. The mesh in the groove area was small to allow the analysis to converge ( Table 2).

Load Application
With the sacrum completely restrained in all directions, a forced displacement of 30 mm to ventral and dorsal direction was applied to the cranial endplate of T8, assuming extension and flexion ( Figure 4b).
An analysis was performed using Patran and MARC (MSC Software, Newport Beach, CA, USA). Sixteen different compression combinations were evaluated, a mechanical analysis of vertebral bodies and discs was performed based on the visual findings, and the highest value for maximum von Mises stress was recorded for each combination. The unit is the megapascal (MPa).

Extension
Maximum stress values at the center of the vertebral body and intervertebral disc from T10 to L2 vertebrae are shown in extension (Table 3 and Figure 5).  Comparisons between the healthy model and DISH flat model showed that stresses were greater on the vertebral body and intervertebral disc from T10 to L2 in extension ( Figure 6) than in the healthy model. Comparisons between the healthy model and DISH flat model showed that stresses were greater on the vertebral body and intervertebral disc from T10 to L2 in extension ( Figure 6) than in the healthy model.  Comparisons between the healthy model and DISH flat model showed that stresses were greater on the vertebral body and intervertebral disc from T10 to L2 in extension ( Figure 6) than in the healthy model.  In the ossified hypertrophic groove 2.5 mm model (D-5/R-2.5 and D-10/R-2.5) (Figure  8), the stress at the center of T12 increased as the groove at the discontinuity became deeper. The maximum stress at the center of each T12 vertebral body is shown in Figures  5 and 8. In extension, the stress on the T11/12 and T12/L1 discs was less than that in the In the ossified hypertrophic groove 2.5 mm model (D-5/R-2.5 and D-10/R-2.5) K (Figure 8), the stress at the center of T12 increased as the groove at the discontinuity became deeper. The maximum stress at the center of each T12 vertebral body is shown in Figures 5 and 8. In extension, the stress on the T11/12 and T12/L1 discs was less than that in the DISH model.  7. In extension analyses of the hypertrophic ossification models T-5 (a) and T-10 (b), the stress on the intervertebral disc at T11-T12 was reduced, whereas stress on the vertebral body was unchanged.
In the ossified hypertrophic groove 2.5 mm model (D-5/R-2.5 and D-10/R-2.5) ( Figure  8), the stress at the center of T12 increased as the groove at the discontinuity became deeper. The maximum stress at the center of each T12 vertebral body is shown in Figures  5 and 8. In extension, the stress on the T11/12 and T12/L1 discs was less than that in the DISH model. In D-5/R-1 and D-10/R-1 (Figure 9), the stress at the center of T12 increased. Maximum stress at the center of each T12 vertebral body is shown in Figures 5 and 10. In the analysis of D-5/R-1, stress on the T12 vertebral body was greater than that of D-5/R-2.5, even with the same depth of ossification. However, with D-10/R-1, the stress on the T12 vertebral body was less than that at D-5/R-1 and at D-10/R-1 in extension. In D-5/R-1 and D-10/R-1 (Figure 9), the stress at the center of T12 increased. Maximum stress at the center of each T12 vertebral body is shown in Figures 5 and 10. In the analysis of D-5/R-1, stress on the T12 vertebral body was greater than that of D-5/R-2.5, even with the same depth of ossification. However, with D-10/R-1, the stress on the T12 vertebral body was less than that at D-5/R-1 and at D-10/R-1 in extension.

Flexion
Maximum stress values at the center of the vertebral body and intervertebral disc from T10 to L2 vertebrae are shown in flexion (Table 4 and Figure 10).

Flexion
Maximum stress values at the center of the vertebral body and intervertebral disc from T10 to L2 vertebrae are shown in flexion (Table 4 and Figure 10).
Comparisons between the healthy model and DISH flat model showed that stresses were greater on the vertebral body and intervertebral disc from T10 to L2 in flexion ( Figure 11) than in the healthy model.  9. In D-5/R-1 (a) and D-10/R-1 (b), the stress at the center of T12 increased at the discontinuity in extension.

Flexion
Maximum stress values at the center of the vertebral body and intervertebral disc from T10 to L2 vertebrae are shown in flexion (Table 4 and Figure 10).    Comparisons between the healthy model and DISH flat model showed that stresses were greater on the vertebral body and intervertebral disc from T10 to L2 in flexion (Figure 11) than in the healthy model. In the T-5 model and T-10 model, the stress on the vertebral body was similar with greater ossification in flexion; no significant differences were observed in stress, except for a decrease in intervertebral disc stress at T11-T12 (Figure 12). In the T-5 model and T-10 model, the stress on the vertebral body was similar with greater ossification in flexion; no significant differences were observed in stress, except for a decrease in intervertebral disc stress at T11-T12 (Figure 12).
(a) (b) Figure 11. Only vertebral bodies and discs are shown for clarity. The right bar is color-coded according to stress (MPa). Stress diagrams of (a) the healthy model and (b) the DISH flat model. Stresses were greater in the flat model in flexion than in the healthy model.
In the T-5 model and T-10 model, the stress on the vertebral body was similar with greater ossification in flexion; no significant differences were observed in stress, except for a decrease in intervertebral disc stress at T11-T12 (Figure 12).
(a) (b) Figure 12. In flexion analyses of T-5 (a) and T-10 (b), the stress on the intervertebral disc at T11-T12 was reduced, whereas that on the vertebral body was unchanged.
In D-5/R-2.5 and D-10/R-2.5 (Figure 13), the stress at the center of T12 increased as the groove at the discontinuity became deeper. The maximum stress at the center of each T12 vertebral body is shown in Figure 13. The results of the flexion analysis were similar to those obtained for the T12 vertebral body. The stress on the T11/12 disc was greater, while that on T12/L1 disc was less than that in the DISH model and the hypertrophic ossification model in flexion. In D-5/R-2.5 and D-10/R-2.5 (Figure 13), the stress at the center of T12 increased as the groove at the discontinuity became deeper. The maximum stress at the center of each T12 vertebral body is shown in Figure 13. The results of the flexion analysis were similar to those obtained for the T12 vertebral body. The stress on the T11/12 disc was greater, while that on T12/L1 disc was less than that in the DISH model and the hypertrophic ossification model in flexion. In D-5/R-1 and D-10/R-1 (Figure 14), the stress at the center of T12 increased. Maximum stress at the center of each T12 vertebral body is shown in Figure 14.
In the analysis of D-5/R-1, stress on the T12 vertebral body was greater than that in D-5/R-2.5, even with the same depth of ossification. However, with D-10/R-1, the stress on the T12 vertebral body was less than that with D-10/R-2.5 and D-5/R-1 in flexion. On the other hand, with D-5/R-1 and D-10/R-1, the stress on the T11/12 and T12/L1 discs increased. The stresses on other vertebral bodies and discs were similar. In D-5/R-1 and D-10/R-1 (Figure 14), the stress at the center of T12 increased. Maximum stress at the center of each T12 vertebral body is shown in Figure 14.
In the analysis of D-5/R-1, stress on the T12 vertebral body was greater than that in D-5/R-2.5, even with the same depth of ossification. However, with D-10/R-1, the stress on the T12 vertebral body was less than that with D-10/R-2.5 and D-5/R-1 in flexion. On the other hand, with D-5/R-1 and D-10/R-1, the stress on the T11/12 and T12/L1 discs increased. The stresses on other vertebral bodies and discs were similar.
In D-5/R-1 and D-10/R-1 (Figure 14), the stress at the center of T12 increased. Maximum stress at the center of each T12 vertebral body is shown in Figure 14.
In the analysis of D-5/R-1, stress on the T12 vertebral body was greater than that in D-5/R-2.5, even with the same depth of ossification. However, with D-10/R-1, the stress on the T12 vertebral body was less than that with D-10/R-2.5 and D-5/R-1 in flexion. On the other hand, with D-5/R-1 and D-10/R-1, the stress on the T11/12 and T12/L1 discs increased. The stresses on other vertebral bodies and discs were similar.

Discussion
DISH is diagnosed based on the continuous calcification of four or more vertebrae, a maintained intervertebral disc height, or no intervertebral joint rigidity [1]. DISH restrains spinal range of motion and causes low back pain [14]. DISH may easily cause fracture in an inflexible spine and spinal damage [15]. This trauma must not be overlooked because the rate of bone union is low and the attenuation of neurological symptoms is difficult to achieve with conservative treatment [2]. However, in contrast to other spine injuries, the diagnosis of DISH is delayed [3,4,16]. One of the factors that make diagnosis difficult is that pre-existent abnormalities on spine radiographs may mask or obscure minimally displaced fractures [3]. Although previous studies reported good outcomes and minimally invasive surgery [17][18][19], the misdiagnosis of fractures and diagnostic delays need to be

Discussion
DISH is diagnosed based on the continuous calcification of four or more vertebrae, a maintained intervertebral disc height, or no intervertebral joint rigidity [1]. DISH restrains spinal range of motion and causes low back pain [14]. DISH may easily cause fracture in an inflexible spine and spinal damage [15]. This trauma must not be overlooked because the rate of bone union is low and the attenuation of neurological symptoms is difficult to achieve with conservative treatment [2]. However, in contrast to other spine injuries, the diagnosis of DISH is delayed [3,4,16]. One of the factors that make diagnosis difficult is that pre-existent abnormalities on spine radiographs may mask or obscure minimally displaced fractures [3]. Although previous studies reported good outcomes and minimally invasive surgery [17][18][19], the misdiagnosis of fractures and diagnostic delays need to be avoided. Because more than 50% of the patients were neurologically intact immediately after the injury, later neurological deterioration occurred in approximately 30% of the total patients [15]. Some previous studies have reported on the risks of trauma of DISH [3,4,16]. However, no report has evaluated the impact of ossification thickness and shape on the spine. Therefore, we considered it important to understand the biomechanics of DISH by analyzing the type of force applied to the vertebral body with different ossification shapes and thicknesses. Although cervical vertebrae are the most frequently injured, we analyzed the thoracolumbar spine transition area because the majority of spinal fractures occur in the thoracolumbar spine [17,20,21].
Recent studies conducted a fracture analysis via 3D-FEM of the vertebral region. Sairyo et al. performed an analysis of lumbar spondylolysis [22]. In addition, Zhao et al. [23] reported the findings of an analysis of compression fractures, while Liang et al. [24] conducted a 3D-FEM analysis of the effects of balloon kyphoplasty; however, only a few intervertebral analyses have been performed to date. Okamoto et al. [25] assessed pressure fractures at different angles of the spine. Nishida et al. [26] and Nakashima et al. [27] conducted a compression fracture analysis of a whole-spine model generated from medical images. However, an analysis and visualization of the stress applied to the vertebral body by changing the ossification pattern around the vertebral body, such as that which occurs in DISH, has not yet been performed. A more detailed understanding of the shapes and discontinuities that are most likely to cause damage may contribute to the prevention of misdiagnoses.
The present study is the first to demonstrate that the stress on the vertebral body was greater with ossification of the anterior longitudinal ligament than with a healthy ligament. The results obtained herein showed the effects of reduced flexibility. Calcification or ossification of the ligament was identified as a contributing factor to increased stress, which may induce vertebral fractures. Okada et al. reported that all of the subjects had anterior and middle column fractures [15]. In this analysis, the greatest increase in stress was from the anterior in both flexion and extension. We found that the partial hypertrophy of ossification did not affect vertebral stress in the ossified area during extension and flexion, where it had an impact on intervertebral discs. The thickness of local ossification did not appear to affect the vertebral body. Furthermore, stresses were more likely to increase at the center of vertebrae with greater ossification and discontinuity such as cases with D-5/R-2.5, D-5/R-1, and D-10/R-2.5. Therefore, if external forces are applied to a case with large ossification and discontinuities, the vertebral groove may fracture. However, this hypothesis was not applied to cases with D-10/R-1. This may be due to the thickness of ossification and the small gap in the curvature, which may have stabilized the vertebral body. On the other hand, the stresses on the cranial and caudal intervertebral discs of discontinuity may have increased during flexion. Okada et al. reported that later neurological deterioration can occur due to a marked instability between the rostral and caudal fused segments [15]. Because some predictions are possible in this analysis, the vertebral body and intervertebral disc around the ossification need to be fully evaluated based on the shape and depth of the bone and the injury mechanism.
The present study had some limitations: Only one type of ossification level, one type of load, one type of bone strength, and one type of spinal alignment were considered. The load was not assumed to be in various alignments, the load of body weight was not analyzed, and the external force was assumed to be a rapid shearing force, not a vertical trauma, such as a buttock injury. Ossification sites are often located to the anterior right of the vertebral body, but not in the middle [28]. Therefore, stresses may change depending on the ossification site. It is not limited to this analysis, but it is a level where increased stress does not necessarily mean fracture. The anterior longitudinal ligament on the caudal side of the vertebral head is ossified and continuous, and posterior elements, such as the vertebral arch, are not turgid. This analysis takes several months to create and analyze a model, which is not at a level that allows for tailor-made diagnosis and treatment. This resolution will require future development of analysis.
However, our analysis showed that the spine is less flexible and more prone to the increased stress of the vertebral body and disc in patients with DISH. Furthermore, the size of ossification and the shape of the groove make it easier to predict levels of increased stress, and an awareness of the possible level of pain and fractures may be useful for clinicians.

Conclusions
The vertebrae of patients with DISH are more susceptible to stress. Furthermore, depending on the morphology of ossification, the stresses on the vertebrae and intervertebral discs differ even with similar loads. Though there is no confirmation concerning effective treatment and diagnosis, an examination of ossification geometry may help surgeons decide the thoracolumbar spine's stress elevated position in patients with DISH, thereby contributing to the understanding of the pathogenesis of pain.  Institutional Review Board Statement: The present study was approved by the Ethics Committee at the Center for Clinical Research, Yamaguchi University Hospital (Ube, Japan; approval No. H28-054). Written informed consent for this study and its publication was obtained from all subjects.
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data presented in this study are available.