Patient-specific finite element analysis for assessing hip fracture risk in aging populations

The femur is one of the most important bone in the human body, as it supports the body’s weight and helps with movement. The aging global population presents a significant challenge, leading to an increasing demand for artificial joints, particularly in knee and hip replacements, which are among the most prevalent surgical procedures worldwide. This study focuses on hip fractures, a common consequence of osteoporotic fractures in the elderly population. To accurately predict individual bone properties and assess fracture risk, patient-specific finite element models (FEM) were developed using CT data from healthy male individuals. The study employed ANSYS 2023 R2 software to estimate fracture loads under simulated single stance loading conditions, considering strain-based failure criteria. The FEM bone models underwent meticulous reconstruction, incorporating geometrical and mechanical properties crucial for fracture risk assessment. Results revealed an underestimation of the ultimate bearing capacity of bones, indicating potential fractures even during routine activities. The study explored variations in bone density, failure loads, and density/load ratios among different specimens, emphasizing the complexity of bone strength determination. Discussion of findings highlighted discrepancies between simulation results and previous studies, suggesting the need for optimization in modelling approaches. The strain-based yield criterion proved accurate in predicting fracture initiation but required adjustments for better load predictions. The study underscores the importance of refining density-elasticity relationships, investigating boundary conditions, and optimizing models through in vitro testing for enhanced clinical applicability in assessing hip fracture risk. In conclusion, this research contributes valuable insights into developing patient-specific FEM bone models for clinical hip fracture risk assessment, emphasizing the need for further refinement and optimization for accurate predictions and enhanced clinical utility.


Introduction
One major challenge with population ageing is the steadily increasing need for artificial joints for older people.Knee and hip replacement operations became some of the most frequently operated surgeries worldwide [1].Their main necessity includes osteoarthritis and, for this study more importantly, osteoporotic fractures [2].These fractures usually result from statically overload or from high forces suddenly applying to the femur and often occur in combination with poor bone mineral density, possibly caused by injuries [3], high age or physical inactivity [1].The mineralization of the bone is strongly correlated to the bone density and directly influencing the strength and load bearing capacity of the joints [4].
The geometrical and mechanical properties must be considered when assessing the risk of bone fracture [5].They strongly depend on bone mineralization and bone shape and vary from person to person, for example due to different age, weight, or sex.To precisely predict the specimen-specific bone properties and their reaction to mechanical loading, there must be a correct geometry reconstruction, an appropriate material assignment and also a realistic replication of the boundary conditions [5].The patient-specific model can then be used in the clinical assessment of bone fracture risk, predicting the fracture load of the femur before a fracture of the bone occurs.This information can be then used for patients who are potentially at risk for osteoporotic fracture.The patients can adapt their motion behaviour to the calculated maximum bearing capacity of the bone and avoid specific movements.This way, fractures can be prevented, and the overall need of replacement surgeries can be reduced.
To assess the hip fracture risk for individual physical activities, from the movement resulting forces on the hip joint need to be estimated.They can be either measured with help of instrumented hip implants or be calculated with gait analysis.In different studies, the resulting force on the femoral head could be determined for important activities like walking, running, and stumbling [2,[6][7][8][9].Peak forces were multiple times higher than the gravitational force of the body weight (BW) and valued about 280%-480% BW during walking, during running up to 615% BW and during stumbling up to 870% BW [10,11].
Previous studies focused on both the numerical analysis with help of the finite element models [5,[12][13][14][15][16][17][18] well as on the experimental in vitro validation of the results obtained by simulations [5,12,14,15,19].Patient-specific finite element models based on CT scans accurately predict the mechanical response and yield load in femurs with metastatic tumors.Additionally, they demonstrate the effectiveness of an elastic quasi-brittle bone description combined with strainbased criteria, which assists in clinical decision-making [20,21].The structural bone models were constructed from quantitative computed tomography (QCT) images.They showed that a linear and isotropic but in homogeneously distributed material is adequate to generate reasonable mechanical models of the human femur [22,23].The two mainly analyzed loading scenarios are the single stance [12-14, 19, 24] and sideways falls [14,15,17,18].For the single stance loading scenario, a static force is applied to the femoral head with an angle of 8 to 20 degree to the femoral shaft in the frontal plane with distal end fixed [13,14].A failure criterion is required to identify the generated models' failure load.It can be distinguished between stress, strain and energy-based failure criteria [25].
This study focusses on the development of patientspecific FEM bone models with possible use in the clinical assessment of hip fracture risk regarding its maximum bearing capacity and the region of crack initiation.For this, ANSYS 2023 R2 was used to estimate the fracture load under a stimulated stance loading condition with a strain-based failure criterion, considering the element-specific strength and stiffness's of the model.

Computed tomography details
For the current study, we considered completely anonymised CT data from healthy four male individual.Table 1 shows the age, height and weight of the individuals considered in this work.Philips Brilliance 64 channel CT scanner was used to obtain the CT DICOM images with slice thickness 0.625 mm.

Femur modelling
The creation of three-dimensional femur geometry is vital in understanding the outcomes of the forces acting on the bone.Figure 1 summarises the workflow for generating femur bones for fracture estimation.(Appendix A shows the patient-specific material assignment from CT scans to 3d model).
Bone geometry was reconstructed from the CT images with the help of 3D Slicer.Manual segmentation is carried out to using a threshold value of 100 HU (Hounsfield Units).The femur model generated from 3D slicer is then exported into CATIA 3D Experience in.STL file format where the bone refinement is carried out, including the smoothening of the bone surface.The improved 3D model of the femur is then exported to the ANSYS 2023 R2 design modeller in.STP file format where the areas are marked for load application.All the eight models trimmed to 120 mm as major area of interest is know the fracture at the neck region of the femur.The geometrical models were meshed with a 3D tetra mesh using hyper mesh with an element size of 2 mm based on the mesh independence study, as shown in figure 2(a).It can be noted that the stress levels below 2 mm show negligible variations, and hence, a mesh with an element size of 2 mm is considered in complete analysis.The number of nodes and elements were given in figure 2(b).In this work, using Bonemat v3.2 software, bone material was considered inhomogeneous and isotropic, relying on the results of [22,23] proving that isotropic material provides sufficient accuracy for strain and stress analyses.The material properties of each element were assigned based on the grey scale of each voxel expressed in the Hounsfield Unit number (HU).The open-source software Bonemat v3.2 and the following empirical equations (1)-(3) were used [17,18,[26][27][28].
is the quantitative computed tomography density, Ash r the ash density, and E refers to the modulus of elasticity.

Loads and boundary conditions
The region for load application was created with the help of a 115 mm offset plane to the distal end so that the height of the application area is in the center of the femoral head and has a constant height of 5 mm for all eight specimens, as shown in figure 3. Referring to [14,15,24], the current study focuses on the simulation of the single stance loading configuration, which is why a static force was applied on the surface of the femoral head at 8°to the diaphysis, with the distal end fully constrained as illustrated in figure 4(b).
Figure 4 shows the young's modulus distribution for data set of 4 R.This has a minimum of 244.9 MPa and a maximum of 8464.7 MPa.The remaining data sets of Young's modulus is given in appendix B.

Results
Finite element analyses were conducted in ANSYS, considering the 'Maximum principal elastic strain' and 'Equivalent stress (von Mises)' as solutions.The failure load for each specimen was determined with the help of the strain yield criterion of 0.61% strain and the maximum stresses within the model were calculated.Table 2 shows the values for the maximum measured bone density, the calculated failure load, and the maximum stress while the failure load was applied to the bone.For reference of the geometry, each bone's volume is also given.The maximum bearing capacity of the bone was determined with the help of the critical strain yield criterion implemented by [14].Failure of the bone is expected as soon as one element exceeds the critical principal strain of 0.61%.To determine the fracture load, a linear strain-force curve was calibrated with two points of arbitrary loads, and then the failure load was calculated by evaluating the curve at 0.61% strain [14].
The calculated failure loads lie between 2.02 and 3.10 kN, what equals to loads 254% to 416% BW of the patients.These values suggest an underestimation of the ultimate bearing capacity of the bones, as the patients did not exhibit any fractures even though the calculated failure loads are considerably below 480% BW and therefore indicate that all the patients could potentially experience fractures during activities like running.
The values of the maximum bone density slightly differ with maximum of 10% deviation, but they are very similar for left and right femur of each patient with maximum of 1% deviation.The ratio between maximum density and failure load is not constant for the patients and varies from 0.41 to 0.62.These values show that bone strength depends on both the maximum bone density and the geometry, so only the mineral density is not a reliable criterion for evaluating bone strength.This observation can also be verified by comparing patients 1 and 2.Even though the maximum density is very similar, the failure loads differ significantly.
The crack initiation region, marked by elements reaching the critical strain value of 0.61%, was in a similar femoral neck region for all specimens.The approximate position is shown in figure 5 for the left femur of patient 3. Stress values of the elements differ significantly between the patients and lie between 31.27 and 51.88 MPa but are rather equal for left and right femur of each patient.In contrast to elements with maximum strain, the corresponding elements were located at the femoral shaft and were loaded with pressure, as shown in figure 6.For patients 1, 3 and 4, only a small difference for the failure load for left and right femur exists.In contrast, for patient 2, a recognizable difference in the properties of the left and right femur can be noted.The failure load deviates 0.63 kN, equal to 24% deviation.The geometrical appearance of the bone cannot easily explain this difference.
It must result from the bone's concrete geometrical shape and density distribution within the bone.The consolidated values of maximum density, failure load and density/load are shown in the graph in figure 7. The maximum density varies across specimens.Specimen 4 R has the highest density, while specimens 3 L and 3 R have relatively lower densities.The failure load is the amount of force applied to the specimen at the point of failure.Specimen 1 L has the highest failure load (3.104 kN), while specimen 2 L has the lowest (2.016 kN).The density/load ratio indicates how efficiently the material can withstand the applied load.Specimen 2 L has the highest ratio (0.62), suggesting it achieves a relatively high load-bearing capacity per unit of density.Table 2. Bone volume, max.Density, failure load, the ratio between density and maximum load, and stress in the yielding element for the eight femurs analysed.

Discussion
As already described in the results, the fracture loads are drastically underestimated by the FE models of the current study.On one side, this can be seen in the values of the fracture loads of 254% to 416% BW of the patients, suggesting that a hip fracture at the femoral neck would occur even during normal walking.On the other side, this hypothesis can be supported by the fracture loads from other studies.The applied loading conditions were identical or similar, but the determined fracture loads were multiple times higher.Cadaveric tests by [24] measured failure loads of 6.3 kN to 11.4 kN equal to 882% to 1413% BW for male patients, respectively.Previously published literature shows cadaveric studies and simulations for two specimens [12].They applied similar boundary conditions but used a principal stress yield criterion and found significantly higher fracture loads of 6020 N and 7120 N for the experiments and 6572 N and 8178 N for the numerical models.The values indicate a slight overestimation of the mechanical properties of the numerical models.Nevertheless, the presented load values are in better accordance with the estimated loads during walking, running and stumbling.
In [15], a very good agreement between experiment and simulation was reached.They measured fracture loads from 6.0 kN to 12.1 kN, while the simulated fracture loads ranged from 4.4 kN to 10.5 kN.
The highest stresses occurred at the femoral shaft in accordance with the stress patterns found by Kheirollahi [25] who applied a force of 250% BW in the direction of the femoral shaft axis.For the single stance, the calculated stresses were higher in the subtrochanteric region than at the femoral neck and the intertrochanteric region of the femur.In addition, they found the highest fracture risk index at the femoral neck which is in line with the in this study determined location of crack initiation.Similarly, in [15], fractures initiated at the femoral neck or the head-neck junction due to tensile strains above 0.73% for all examined specimens.In contrast to this, the results obtained by Laz et al [13] showed discrepancies.The force was applied at a very similar surface area but at 20°to the femoral shaft axis.The highest stresses occurred in the region of the femoral neck, where also the highest risk of fracture was supposed.The simulated stress-induced failure loads reached from 15.2 to 25.2kN.They exceeded the results of other studies by far, which makes the accuracy of the estimated failure loads and used stress failure criterion questionable.As also in clinical cases most fractures are femoral neck and intertrochanteric hip fractures and only the uncommon subtrochanteric fracture is expected to be stress-induced [15,25,29], the maximum strain yield criterion as implemented in the current study seems to be more accurate than a stress failure limit.
As already documented by [14], this study provided a good estimation of the fracture pattern but does not generate accurate values for the failure loads.They applied the same boundary conditions and used the maximum strain yield criterion of 0.61%.The fracture location was also at the femoral neck and matched closely between experiment and simulation.However, the simulated fracture loads for the stance loading condition were approximately 50% smaller than those measured during the experiments and reached between 1954N and 3927 N, comparable to the current results.This underestimation of the failure loads can be reasoned with applying the strain yield criterion to one element, leading to an increased sensitivity of the surface elements [14].However, as the fracture region was correct, the strain yield criterion seems suitable but could be increased to a larger threshold value as done by other studies [15,30].Composite bone plates with varied fibers and orientations can be used to reduce stress shielding in femur fracture fixation, aiming for improved biomechanical performance compared to conventional metal plates [31].The incidence of postoperative periprosthetic femoral fractures necessitates tailored management strategies according to femoral stem characteristics.Future research should prioritize identifying implantrelated risk factors and optimizing concurrent metabolic bone disease treatment [32].
A prior optimization of the models is required to assess the risk of hip fracture in clinical cases.Several factors influence the simulation and need to be investigated.This is, very importantly, the used densityelasticity relationship that strongly interacts with the obtained results [5].To do so, a combined study of in vitro testing and simulation, as done by [5], is required to assign the material properties more accurately to the models.Moreover, the boundary conditions provided in the human body must be investigated and adequately reproduced to correlate the calculated fracture loads with the in vivo measured or calculated ultimate strength of the bone.Lastly, the strain-based yield criterion provides a correct position for crack initiation, but its value needs to be adapted to predict better the failure loads.

Limitations
The study only considered CT data from four healthy male individuals.This limited sample size may not capture the full variability in bone geometry and mechanical properties among different populations, such as females or individuals with osteoporosis.The study makes certain assumptions and simplifications in the modeling process, such as ignoring soft tissue effects and considering only one loading scenario (single stance).These assumptions may not fully capture the complexity of hip fracture mechanisms in real-world scenarios [33].While the study aims to assess hip fracture risk, the clinical relevance of the predicted fracture loads to real-world fracture occurrence and prevention strategies is not fully explored.Further research is needed to translate the findings into practical clinical applications for fracture prevention and management.Addressing these limitations through larger sample sizes, more realistic material models, improved boundary conditions, direct validation against experimental data, exploration of different failure criteria, and consideration of clinical implications can enhance the robustness and applicability of the current research findings.

Conclusion
In conclusion, this study delves into the critical realm of hip fracture risk assessment through the development of patient-specific finite element models (FEM) for femur bones.The research, focused on a cohort of healthy male individuals, utilized CT data to construct accurate geometric and mechanical representations, employing ANSYS 2023 R2 software to simulate fracture loads under single stance loading conditions.The findings reveal a significant underestimation of the ultimate bearing capacity of bones, hinting at potential fractures even during routine activities.Notably, discrepancies between simulation results and prior studies underscore the necessity for optimization in modeling approaches.The study emphasizes the complexity of bone strength determination, exploring variations in bone density, failure loads, and density/load ratios among specimens.The strain-based yield criterion proves effective in predicting fracture initiation but requires adjustments for more precise load predictions.The research underscores the importance of refining density-elasticity relationships, investigating boundary conditions, and optimizing models through in vitro testing for enhanced clinical applicability in assessing hip fracture risk.As the aging global population amplifies the demand for artificial joints, especially in knee and hip replacements, the insights gained from this research are invaluable.Despite the challenges and discrepancies encountered, the study lays the foundation for future refinements and optimizations in modeling approaches, with the ultimate goal of providing more accurate predictions for enhanced clinical utility.Further interdisciplinary research combining in vitro testing and simulation, along with a comprehensive exploration of density-elasticity relationships, is crucial to advancing our understanding and prediction of hip fracture risks in clinical scenarios.

Figure 1 .
Figure 1.Workflow for fracture analysis of bone.

Figure 4 .
Figure 4. (a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 4 R.

Figure 3 .
Figure 3. (a) Femur model with marked region for load application.(b) Applied boundary conditions for single stance loading case.

Figure 7 .
Figure 7. Variations of Failure load, density and density/load for all the specimens.

Figure
Figure B1.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 1LR.

Figure B2 .
Figure B2.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 1 R.

Figure B4 .
Figure B4.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 2 R.

Figure B5 .
Figure B5.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 3 L.

Figure B6 .
Figure B6.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 3 R.

Figure B7 .
Figure B7.(a) Young's model distribution for data set 4 R. (b) Cross section with Young's model distribution for data set 4 L.

Table 1 .
Patient details consider for this work.