An energy approach describes spine equilibrium in adolescent idiopathic scoliosis

Abstract

The adolescent idiopathic scoliosis (AIS) is a 3D deformity of the spine whose origin is unknown and clinical evolution unpredictable. In this work, a mixed theoretical and numerical approach based on energetic considerations is proposed to study the global spine deformations. The introduced mechanical model aims at overcoming the limitations of computational cost and high variability in physical parameters. The model is constituted of rigid vertebral bodies associated with 3D effective stiffness tensors. The spine equilibrium is found using minimization methods of the mechanical total energy which circumvents forces and loading calculation. The values of the model parameters exhibited in the stiffness tensor are retrieved using a combination of clinical images post-processing and inverse algorithms implementation. Energy distribution patterns can then be evaluated at the global spine scale to investigate given time patient-specific features. To verify the reliability of the numerical methods, a simplified model of spine was implemented. The methodology was then applied to a clinical case of AIS (13-year-old girl, Lenke 1A). Comparisons of the numerical spine geometry with clinical data equilibria showed numerical calculations were performed with great accuracy. The patient follow-up allowed us to highlight the energetic role of the apical and junctional zones of the deformed spine, the repercussion of sagittal bending in sacro-illiac junctions and the significant role of torsion with scoliosis aggravation. Tangible comparisons of output measures with clinical pathology knowledge provided a reliable basis for further use of those numerical developments in AIS classification, scoliosis evolution prediction and potentially surgical planning.

Appendices

Appendix 1: Rotation matrix

Rotation matrix $R_i$ :

$$\begin{aligned} \begin{pmatrix} \cos \theta _{i} \cos \alpha _{i} &{} -\sin \theta _{i} \cos \alpha _{i} &{} \sin \alpha _{i} \\ \sin \theta _{i} \cos \varphi _{i} + \cos \theta _{i} \sin \alpha _{i} \sin \varphi _{i} &{} \cos \theta _{i} \cos \varphi _{i} - \sin \theta _{i} \sin \alpha _{i} \sin \varphi _{i} &{} - \cos \alpha _{i} \sin \varphi _{i} \\ \sin \theta _{i} \sin \varphi _{i} - \cos \theta _{i} \sin \alpha _{i} \cos \varphi _{i} &{} \cos \theta _{i} \sin \varphi _{i} + \sin \theta _{i} \sin \alpha _{i} \cos \varphi _{i} &{} \cos \alpha _{i} \cos \varphi _{i} \end{pmatrix} \end{aligned}$$

(17)

Appendix 2: Adjoint method

The derivative operator $\nabla _{\mathbf {p}}$ used on any function $h({\mathbf {u}}_\text {eq}({\mathbf {p}}), {\mathbf {p}})$ is defined as:

$$\begin{aligned} \nabla _{\mathbf {p}} h = \frac{\hbox {d}h}{\hbox {d}{\mathbf {p}}} \end{aligned}$$

(18)

Therefore,

$$\begin{aligned} \nabla _{{\mathbf {p}}} f = \frac{\partial f}{\partial {\mathbf {u}}_\text {eq}} \frac{\hbox {d} {\mathbf {u}}_\text {eq}}{\hbox {d} {\mathbf {p}}} + \frac{\partial f}{\partial {\mathbf {p}}}. \end{aligned}$$

(19)

Calling ${\mathbf {g}}({\mathbf {u}}, {\mathbf {p}}) = (\nabla _{\mathbf {u}} V)({\mathbf {u}}, {\mathbf {p}})$, the definition of ${\mathbf {u}}_\text {eq}$ for every ${\mathbf {p}}$ gives:

$$\begin{aligned} {\mathbf {g}}({\mathbf {u}}_\text {eq}({\mathbf {p}}), {\mathbf {p}})=0 \quad \text {therefore} \quad \nabla _{\mathbf {p}}({\mathbf {g}}({\mathbf {u}}_\text {eq}({\mathbf {p}}), {\mathbf {p}})) = 0, \end{aligned}$$

(20)

leading to the following equation:

$$\begin{aligned} \frac{\hbox {d} {\mathbf {u}}_\text {eq}}{\hbox {d} {\mathbf {p}}} = - \left( \frac{\partial {\mathbf {g}}}{\partial {\mathbf {u}}_\text {eq}}\right) ^{-1} \frac{\partial {\mathbf {g}}}{\partial {\mathbf {p}}}. \end{aligned}$$

(21)

Using this result in Eq. (19) gives:

$$\begin{aligned} \nabla _{{\mathbf {p}}} f = - \frac{\partial f}{\partial {\mathbf {u}}_\text {eq}} \left( \frac{\partial {\mathbf {g}}}{\partial {\mathbf {u}}_\text {eq}}\right) ^{-1} \frac{\partial {\mathbf {g}}}{\partial {\mathbf {p}}} + \frac{\partial f}{\partial {\mathbf {p}}} \end{aligned}$$

(22)

The vector $- \frac{\partial f}{\partial {\mathbf {u}}_\text {eq}} \left( \frac{\partial {\mathbf {g}}}{\partial {\mathbf {u}}_\text {eq}}\right) ^{-1}$, often called $\varvec{\lambda }$, is the adjoint vector. It is the solution of a linear system, faster to solve than an explicit finite difference calculation to access the gradient $\nabla _{\mathbf {p}}\, f$.

Appendix 3: Sensitivity study

In the proposed methodology, the spine balance was found by minimizing the total mechanical energy using EOS^® medical images from patient follow-up as input data. We assessed the impact of uncertainties of clinical data numerization on numerical prediction.

The measurements errors were evaluated using ten numerizations of a single clinical image. The locations and orientations of seventeen vertebral bodies have been computed for each numerization. Inspired by Bayesian methodology, the geometry space was described with a probability law chosen to be normal, i.e., defined by mean value and standard deviation. The first and last decile of the distribution provided the envelops of the clinical geometry. The maximum distance between the envelops was evaluated at 5% of the maximum displacement from vertical spine in both frontal and sagittal direction.

The computation cost for the uncertainties propagation through the inverse problem was prohibitive. Therefore, we investigated the propagation of hypothetical uncertainties on the parameters, through the direct problem. The parameters uncertainties were chosen to be independent and following a normal distribution with a fixed mean and an arbitrary initial standard deviation (few percents of the mean). The deciles distribution of the computed equilibrium geometry ${\mathbf {u}}_\text {eq}$ was then compared with the deciles of clinical data from image numerization. The parameters standard deviation was iteratively updated to obtain a good match between clinical measurements envelops and equilibrium geometry envelops. After computation, the discrepancies between clinical envelops and equilibrium envelops was lower than 6% and was obtained with uncertainties on parameters characterized by standard deviation of 5% of the mean values.