Advanced multi-scale mechanobiological simulations enable the distinction between healthy and pathological bone growth patterns

Abstract

Introduction & Purpose

Bones adapt its shape and strength in response to mechanical loading. Tissue-level mechanical stresses, i.e., hydrostatic compressive stress and octahedral shear stress, regulate biological events on the cellular and molecular level, and thereby impact tissue histomorphology and changes in bony geometry (Carter & Beaupré, 2007). Mechanobiological multi-scale simulations based on gait analysis data, musculoskeletal (MSK) modelling and finite element (FE) analysis enable to estimate tissue-level bone loads and predict growth patterns (Shefelbine & Carter, 2004). Such simulations have the potential to inform clinical decision-making in patients with pathological bone deformations, e.g. increased femoral anteversion angle (AVA), which are common in children with cerebral palsy (CP). However, the time-consuming process of creating subject-specific models has previously restricted studies to small sample sizes (n = 1-4 participants), thereby limiting the generalizability of findings (Carriero et al., 2011; Yadav et al., 2021). Moreover, the simulations have not been validated directly due to a lack of longitudinal experimental data, which hinders the widespread use and implementation in clinical settings. Over the last years, we focused our research efforts to overcome these hurdles. We developed a workflow to streamline the process to create subject-specific FE model (Koller et al., 2023). In this study we aimed to use this workflow and 1) analyze growth patterns in a comparable large sample of typically developing (TD) children and children with CP and 2) validate the growth predictions, i.e. compare predicted with measured bone growth.

Methods

Magnetic resonance images (MRI) and three-dimensional motion capture data including marker trajectories and ground reaction forces of 13 TD children (10 ± 2.2 years old, mass: 36.8 ± 9.5 kg) and 12 children with CP (10.4 ± 3.8 years old, mass: 30.1 ± 10.8 kg) were analyzed for this study. All participants walked without walking aids and with a self-selected speed. Details of the data collection are described in previous publications (Kainz et al., 2017; Koller et al., 2023). From 10 TD children we collected the MRI and motion capture data twice with approximately two years between occasions.

MRI images were used to create subject-specific MSK models (Veerkamp et al., 2021) based on an OpenSim model (Kaneda et al., 2023), which allows to estimate medial and lateral knee joint contact forces (JCF) as well as the patellofemoral JCF. Each participant’s model and the corresponding gait analysis data were used to calculate joint angles, muscle forces and JCFs using inverse kinematics, static optimization and joint reaction load analyses with OpenSim 4.2, respectively. For each participant, muscle forces and JCF of two representative steps (left and right) were chosen as loading condition for the FE simulations and each side was analyzed separately.

Each femur was segmented using 3D Slicer and the geometry was used to calculate the bending of the femoral shaft. Our previously developed tool was used to create subject-specific hexahedral meshes with several layers of elements aligned with the growth plates (Koller et al., 2023). Nine load instances were selected based on the JCF peaks and the valley in-between during the stance phase. JCFs and muscle forces (n = 26) acting on the femur at these timepoints were used as loading conditions for FE analysis. Linear elastic material properties were assigned to the different parts of the femur and FEBio3 was used to calculate principal stresses within the growth plates.

The growth rate due to mechanical loading was estimated as the osteogenic index (OI), which was calculated based on the hydrostatic compressive stress and octahedral shear stress for each element within the growth plates (Stevens et al., 1999; Yadav et al., 2021). Positive and negative OI values indicate regions where growth is likely to be promoted or inhibited, respectively. Image comparison methods were employed to identify variability between OIs (Bradski, 2000).

The knee flexion angle alters the orientation of knee JCFs in respect to the distal growth plate resulting in different induced stress regimes. Hence, we divided the children with CP in two groups: one with normal knee flexion angles (n = 15 femurs, CP_normal) and one with increased knee flexion (n = 8 femurs, CP_{high_knee_flexion}) during the stance phase. For aim 1, we compared the femoral geometry and OI between CP groups and the TD children.

For aim 2 (validation), growth predictions were performed for all femurs where data of two occasions was available (n = 20 femurs). Monte-Carlo-analysis (n = 1,320 for each femur) were performed to estimate which mechanobiological model parameters lead to the most realistic simulations. Linear regression analyses were used to identify whether measured development of AVA can be predicted by the mechanobiological multi-scale simulations.

Results

The bending of the femoral shaft was significantly different between all groups (Figure 1B). The representative reference OI distribution from the TD femurs showed a ring-shape at the proximal growth plate. Within the CP cohort, image comparison revealed higher inter-subject variability (p < 0.001) whereas some participants show ring-shape distribution similar to the TD cohort and others linear-gradient shapes with high values on the lateral side of the growth plate. At the distal growth plate, OI distribution for TD femurs showed highest values in the posterior notch region and on the medial-anterior edge. The OI distribution of the CP_normal group had its peak values in similar areas as the TD group but with additional peaks at the anterior-lateral edge. The OI distribution generated from the CP_{high_knee_flexion} group showed a linear gradient from high to low values from anterior to posterior side (Figure 1D).

Between data collection sessions, TD children grew 13.5 ± 3.3 cm and gained 8.2 ± 3.2 kg of body mass during the two years. Participants’ AVA changed between -13.1° and 11.8° (mean: -1.3 ± 5.8°) between sessions. Analysis of kinematics, muscle and JCFs did not show differences between those experiencing an increase of AVA to those that showed a decrease. Multi-scale predictions and measurements of AVA development showed significant linear correlation (p < 0.001) with an explanatory power of R² = 0.5 (Figure 1C).

Discussion

At the distal growth plate, the OI distribution showed more growth in the anterior compared to the posterior region in the CP_{high_knee_flexion} compared to the CP_normal and TD groups. Progressive promoted growth in the anterior compartment over a longer period could lead to higher anterior bending of the femoral shaft. Indeed, our geometrical analysis revealed a more bended femoral shaft in the CP_{high_knee_flexion} group compared to other groups.

The variability of the OI at the proximal growth plate within the CP cohort was higher compared to the TD group indicating that some CP participants are likely to experience normal growth whereas in others, the proximal femur will develop pathologically. This agrees with clinical observation where some children with CP develop deformities while others don’t.

Despite the fact that no significant differences were found in joint kinematics, femoral loading and growth plate orientation between TD children with different growth patterns, multi-scale simulations were sensitive enough to identify differences and predict AVA development with reasonable accuracy.

Conclusion

Femoral growth is influenced by a complex interplay between gait pattern, femoral morphology and internal loading on a tissue-level. Our results showed that multi-scale simulations are able to discriminate between different growth patterns and predict growth trends in agreement with experimental observations.  In order to increase our confidence in the simulations and pave the road to in-silico informed clinical decision-making in the near future, longitudinal simulation studies including a larger sample size and individuals with pathological growth are needed to be conducted.

References

Bradski, G. (2000). The OpenCV library. Dr. Dobb’s Journal of Software Tools.

Carriero, A., Jonkers, I., & Shefelbine, S. J. (2011). Mechanobiological prediction of proximal femoral deformities in children with cerebral palsy. Computer Methods in Biomechanics and Biomedical Engineering, 14(3), 253–262. https://doi.org/10.1080/10255841003682505

Carter, D. R., & Beaupré, G. S. (2007). Skeletal function and form: Mechanobiology of skeletal development, aging, and regeneration. Cambridge University Press.

Kainz, H., Hoang, H. X., Stockton, C., Boyd, R. R., Lloyd, D. G., & Carty, C. P. (2017). Accuracy and reliability of marker-based approaches to scale the pelvis, thigh, and shank segments in musculoskeletal models. Journal of Applied Biomechanics, 33(5), 354–360. https://doi.org/10.1123/jab.2016-0282

Kaneda, J. M., Seagers, K. A., Uhlrich, S. D., Kolesar, J. A., Thomas, K. A., & Delp, S. L. (2023). Can static optimization detect changes in peak medial knee contact forces induced by gait modifications? Journal of Biomechanics, 152, Article 111569. https://doi.org/10.1016/j.jbiomech.2023.111569

Koller, W., Gonçalves, B., Baca, A., & Kainz, H. (2023). Intra- and inter-subject variability of femoral growth plate stresses in typically developing children and children with cerebral palsy. Frontiers in Bioengineering and Biotechnology, 11, Article 1140527. https://doi.org/10.3389/fbioe.2023.1140527

Shefelbine, S. J., & Carter, D. R. (2004). Mechanobiological predictions of femoral anteversion in cerebral palsy. Annals of Biomedical Engineering, 32(2), 297–305. https://doi.org/10.1023/B:ABME.0000012750.73170.ba

Stevens, S. S., Beaupré, G. S., & Carter, D. R. (1999). Computer model of endochondral growth and ossification in long bones: Biological and mechanobiological influences. Journal of Orthopaedic Research, 17(5), 646–653. https://doi.org/10.1002/jor.1100170505

Veerkamp, K., Kainz, H., Killen, B. A., Jónasdóttir, H., & van der Krogt, M. M. (2021). Torsion Tool: An automated tool for personalising femoral and tibial geometries in OpenSim musculoskeletal models. Journal of Biomechanics, 125, Article 110589. https://doi.org/10.1016/j.jbiomech.2021.110589

Yadav, P., Fernández, M. P., & Gutierrez-Farewik, E. M. (2021). Influence of loading direction due to physical activity on proximal femoral growth tendency. Medical Engineering & Physics, 90, 83-91. https://doi.org/10.1016/j.medengphy.2021.02.008