Separate and Combined Effects of Geometrical and Mechanical Properties Changes Due to Aging on the Femoral Strength in Men and Women

Aging, from 40 to +80 years old, causes geometrical and mechanical properties changes in the proximal femur. The subperiosteal width expands faster in men compared to women during aging, while the cortical thickness varies unequally in each sector and differently between men and women. Another change which occurs during aging is bone mechanical properties such as stiffness and ultimate strains. Numerical analysis allows us to study the potential effects of each of the age-dependent changes on the fracture forces separately and combined. We investigated the effects of the geometrical and bone mechanical properties changes due to aging on the femoral strength during a common falling scenario using a transverse isotropic continuum damage model. First, the femur model was adapted from a previously developed human body model named THUMS v4.02. Then, three sets of models were developed to address each of the changes separately and combined for both sexes. We found that the fracture forces in women are on average 1500 N less than in men of the same age. The age-dependent geometrical changes increased the fracture forces in men (25 N/decade), whereas it reduced the fracture forces by 116 N/decade in women. The mechanical properties changes reduced the fracture forces in men more than in women (354.5 N/ decade vs. 225.4 N/decade). When accounting for both geometrical and mechanical properties changes due to aging, the fracture forces decreased by 10.7% of the baseline in women per decade compared to 7.2% per decade in men.


INTRODUCTION
Adults older than 65 are at risk of a fall, and this risk increases as they get older (Ambrose et al., 2013). During 2014, 25% of older adults in the United States experienced falling at least once (Bergen et al., 2016). Falls are the leading cause of accidental death and nonfatal injuries in this age group . Nonfatal and fatal falls cost $50 billion, and falling is the costliest nonfatal injury in the US (Florence et al., 2018). Houry et al. (2016) estimated growth in the population of older adults, which without preventive plans, could lead to a 100% increase in the number of injuries caused by a fall in the year 2030.
Falls, commonly to the sideways, cause more than 95% of all hip fractures (Hayes et al., 1993;Parkkari et al., 1999;Nasiri Sarvi and Luo 2017). Hip fractures lead to hospitalization and cause long-term consequences on mobility and independence of the patients (Florence et al., 2018). Elderly women are at a higher risk of hip fracture comparing to men at a similar age (Courtney et al., 1995;Silva 2007;Ito et al., 2011), and the fracture risk increases exponentially as men and women get older (Beck et al., 2000;Ito et al., 2011).
The geometry of the bone, bone mechanical properties, and loading mode are three factors that affect the femoral strength (Silva 2007). The first two factors are changing with aging, and the loading mode represents an extrinsic factor. The geometrical parameters of the bone, such as the cortical thickness, cortical cross-sectional area, and cross-sectional moment of inertia, are negatively correlated with aging (Ito et al., 2011). In other words, the endosteal resorption and periosteal expansion cause geometrical changes of the bone. In general, the cross-sectional area of the femoral neck is larger in men comparing to women; and these differences increase as they get older (Beck et al., 2000;Carpenter et al., 2011). Several studies suggested that the cortical thicknesses are changing dissimilarly in each quadrant of the femoral neck (Mayhew et al., 2005;Poole et al., 2010;Johannesdottir et al., 2011;Poole et al., 2012;Carballido-Gamio et al., 2013;Johannesdottir et al., 2013). A longitudinal study by Johannesdottir et al. (2013) suggested that the average superior cortical thickness gets thinner at a higher rate than the inferior thickness in a 5-year period for both sexes. Moreover, the thickness decreases slower in men comparing to women. Similar behavior was reported in several cross-sectional studies (Poole et al., 2010;Ito et al., 2011;Johannesdottir et al., 2011;Khoo et al., 2019). Mayhew et al. (2005) have reported the corresponding cortical thicknesses in 16 sectors of the femoral neck. The average superior cortical thickness decreases with aging in men and women (Poole et al., 2010;Ito et al., 2011;Johannesdottir et al., 2011;Khoo et al., 2019), while the average inferior cortical thickness remains constant with aging. Mechanical properties of the cortical and trabecular bone are changing with age, and these changes are site-dependent (Osterhoff et al., 2016). Unlike as observed for vertebrae, trabecular bone in the femur merely has a minor contribution to the biomechanical strength of the whole-bone compared to the cortical bone (Holzer et al., 2009;Osterhoff et al., 2016). In the femoral neck, the cortical bone loses stiffness in terms of Young's and shear modulus by 3 and 4% per decade of age, respectively. Ultimate strains decrease about 5-10% of the initial value per decade while yield properties of the bone do not change significantly with age (Reilly et al., 1974;Reilly and Burstein 1975;Martin 1993;Dokko et al., 2009;Kutz 2009;Mirzaali et al., 2016;Osterhoff et al., 2016). Those changes, combined, cause a transition from a ductile bone to a more brittle one with age (Ott et al., 2010). Trabecular bone is experiencing similar mechanical properties changes in the femoral neck where the modulus is decreasing approximately 5% per decade (Mosekilde et al., 1987;Kutz, 2009). In contrary to the geometrical changes, mechanical properties changes of bone are sexindependent (Martin 1993;Morgan and Bouxsein 2008;Osterhoff et al., 2016).
Several whole-bone experiments and finite element models are designed to investigate the aging effect on the femoral strength. Age, sex, bone mineral density (BMD), loading rate, loading configuration, and individual geometrical characteristics are widely used to explain the variation of the femoral strength in test subjects (Silva 2007;Morgan and Bouxsein 2008;Dragomir-Daescu et al., 2011;Grassi et al., 2012;Panyasantisuk et al., 2018). A recent study by Dragomir-Daescu et al. (2018) has concluded that the aBMD, sex, and age are sufficient parameters to clinically evaluate the femoral strength. In addition to them, they also explored the effect of loading rate and neck-shaft angle. They discovered that the inclusion of the two latter variables would not considerably improve the prediction of the femoral strength (Dragomir-Daescu et al., 2018). Several other clinical and whole-bone experiments indicate the effect of bone size and mechanical properties of the bone on the femoral strength (Bouxsein 2006;Morgan and Bouxsein 2008).
The majority of current computational studies on the strength of femur are focused on developing subject-specific finite element models, which are typically using isotropic mechanical properties to evaluate the femur strength (Dall'Ara et al., 2013;Rezaei et al., 2017). While those models are capable of predicting the subject-specific femoral strength accurately, they are not developed to predict the effect of changes that occur to the bone of average men and women due to aging. There has been no work published so far to the best of authors' knowledge, which explicitly studies the separate and combined effects of the geometrical and mechanical properties changes due to aging in each of the sexes. Moreover, many studies (Reilly et al., 1974;Reilly and Burstein 1975;Ashman et al., 1984;Cezayirlioglu et al., 1985;Rho et al., 1995;Derikx et al., 2011;Mirzaali et al., 2016;Falcinelli and Whyne 2020) have shown that the yield and ultimate strains of the cortical bone are asymmetric which in turn could affect the femoral strength predictions and fracture initiation. The subject-specific models are also getting computationally expensive to the extent that it becomes impractical to implement them into whole-body models. In the current study, we aim to investigate the separate and combined effect of geometrical and mechanical properties changes on fall-induced femur fractures using a transverse isotropic continuum damage model. Simultaneously, it is intended to keep the simplicity of the model while improving the accuracy of bone modeling to allow direct application of this method in the whole-body models. We applied the geometrical and mechanical properties change for 4 decades of age, from 40 to 80 years old with the aim to evaluate the effect of each of the age-dependent changes on the femoral strength.

MATERIALS AND METHODS
Three distinct sets of numerical models were developed to investigate the separated and combined effects of key changes in the geometry and mechanical properties on the strength of femur. The first set accounted only for age-dependent geometrical changes (GeomAge), the second set accounted only for age-dependent mechanical properties changes (MechAge), and the third set accounted for both changes in 4 decades of age (MechGeomAge), from 40 to 80 years old.

Baseline Model
The geometry of the baseline model for 40 years old men and women was obtained from the total human model for safety (THUMS) v4.02. It represents the 50th percentile adult man with the height and weight of 77 and 1.75 m, respectively (Iwamoto et al., 2002). Albeit, the initial model needed several modifications to represent elderly men and women. Foremost, the cortical thickness was uniform around the femoral neck in the initial model, whereas the cortical thickness naturally varies in each sector of the femoral neck (Mayhew et al., 2005;Johannesdottir et al., 2013;Khoo et al., 2019). Additionally, the subperiosteal width of the narrowest point of the neck is different among men and women (Beck et al., 2000;Duan et al., 2003), which indicated the significance of developing two distinctive baseline models for men and women.
The geometry of the initial THUMS model is based on a healthy man, and the subperiosteal width is in close agreement with the reported value from the national health and nutrition examination survey (NHANES III) database (Beck et al., 2000;Duan et al., 2003). On the other hand, the baseline model for women required scaling down to the representative subperiosteal width following the NHANES III (Beck et al., 2000). The homogenous scaling factor for this adjustment was around 0.859. Afterward, the baseline model for men and women was adjusted for the corresponding cortical thickness in each sector of the femoral neck. The cortical thicknesses in each sector were calculated using the combination of the average cortical thicknesses for each age group from the NHANES III database (Beck et al., 2000;Duan et al., 2003) and the relative thicknesses in each sector of the femoral neck from the study by Mayhew et al. (2005). Since the average inferior cortical thicknesses in the study by Mayhew et al. (2005) was not in agreement with similar studies (see Supplementary Material), we decided to transform the reported cortical thicknesses of each sector (d sector i ) to a relative thickness by dividing those thicknesses by the maximum reported thickness (D Max reported ). Then, we used the average of the relative cortical thicknesses and the NHANES III average cortical thicknesses (D avg, NHANES III ) for each age group to evaluate the new maximum cortical thickness (D NewMax ) (Eq. 1): (1) Finally, the relative cortical thicknesses of each sector were transformed back to absolute thicknesses using the new maximum thickness (Beck et al., 2000;Mayhew et al., 2005;Figures 1A,B).
The mechanical properties of both cortical and trabecular bones were modified. The material model of the cortical bone was altered to a transverse isotropic continuum damage model since it better can represent the cortical bone than an isotropic model (Reilly et al., 1974;Reilly and Burstein 1975;Dong and Guo 2004;de Bakker et al., 2009) ; and it accounts for the asymmetric yield and ultimate strains (Reilly et al., 1974;Reilly and Burstein 1975;Ashman et al., 1984;Cezayirlioglu et al., 1985;Rho et al., 1995;Derikx et al., 2011;Falcinelli and Whyne 2020;Mirzaali et al., 2016) in tension, compression, and shear in respective directions. Mechanical properties of the cortical bone were adapted from the experiments by Reilly et al. (1974) and Reilly and Burstein (1975 ,  Table 1). Those properties were assumed to be similar for both men and women.
The transverse isotropic elastic model was combined with a linear damage model. Six damage parameters were applied independently to the stiffness parameters in respective directions. The damage initiates as an element reaches the yield strain in the respective direction. Subsequently, it linearly decreases the stiffness in the corresponding direction until the ultimate strain is reached. In each timestep, the damage parameters were automatically calculated using Eq. 2: where ( ) + is the positive part of: ; ϵ s k and ϵ c k are yield and ultimate strains in the respective directions. Thexstiffness parameters were set to 1% of their initial value as the strain reaches its ultimate value. It adequately secured the integrity of the mesh. Moreover, it alleviated the demand for a smoother mesh, which might increase the computational costs for the whole-body model. The update to the mechanical properties of the trabecular bone was limited to the stiffness and ultimate strain values. The material model remained isotropic elastic-plastic with continuum damage mechanics. The geometry of the cortical bone was divided into two distinct sections. The shaft and neck region were assumed to be transverse isotropic. The principal axis of the transverse plane is respectively along the shaft and neck hollow ( Figure 1C). Contrary to the shaft and neck, the femoral head typically undergoes multidirectional loading; therefore, it was assumed to have an isotropic continuum damage model with similar mechanical properties to the average value of the shaft and neck regions. The cortical bone was meshed with selectively reduced hexahedron elements, and the trabecular bone was meshed using tetrahedron elements (average dimension of elements 1.7 mm, and the total number of 32,823 nodes and 127,821 elements). A convergence study on the mesh was performed using a one-step division of the cortical bone elements, and remeshing of the trabecular bone according to the new surrounding cortical bone elements. It was found that the maximum forces differed approximately 10% and the refined model was utilized in this study.

Aging Effects
Two geometrical parameters were varied in the GeomAge set, which corresponds to the subperiosteal expansion and the endocortical resorption. Two distinct models were developed for men and women since those changes are different between sexes. First, the subperiosteal width of the baseline model was scaled to the corresponding width for each decade (Beck et al., 2000;Duan et al., 2003). The width was increased by 0.08 and 0.028 mm/year for men and women, respectively (Duan et al., 2003). Afterward, the cortical thickness of each sector of the femoral neck was modified to the corresponding thicknesses obtained from the combination of the NHANES III database Beck et al., (2000) and Mayhew et al. (2005) as explained for the baseline models. The cortical thickness of the femoral neck was assumed to be constant along the neck axis in each sector. Moreover, it changed linearly during aging. The angle in Figures 1A,B corresponds to the rotation in the clockwise direction, starting from the superior surface of the femur.
In the MechAge set, mechanical properties of the bone, such as stiffness, were changed due to aging. The Young's and shear modulus of the cortical bone was reduced by 3 and 4% per decade of age, respectively (Burstein et al., 1976). The ultimate strain was decreased by 10% of the initial value per decade, and yield properties were remained constant during aging (Reilly et al., 1974;Reilly and Burstein 1975;Martin 1993;Dokko et al., 2009;Kutz 2009;Mirzaali et al., 2016). Likewise, the modulus of elasticity and yield stresses of the trabecular bone were reduced 17 and 1.09 MPa, respectively, per decade of age (Mosekilde et al., 1987;Keaveny, 1998). The changes were assumed to be linear and the same for all directions. It is noteworthy that the age in the current study is a representative age which means the bone strength of a 40 years-old osteoporotic subject is substantially comparable with an older representative age group in the current simulations.
The geometry from the GeomAge set was combined with the corresponding mechanical properties from the MechAge set to produce the MechGeomAge set. Since the GeomAge consisted of both sexes, the MechGeomAge could explore the potential differences between men and women as well. The MechGeomAge conveniently represented the general aging effects on the femur for both men and women.

Fall Simulations
The femur was positioned in the sideways fall experimental setup in all three sets. Initially, the femoral shaft was rotated 10°with respect to the horizontal plane, and the femoral neck was rotated internally to 15° (Courtney et al., 1994;Dragomir-Daescu et al., 2011). Next, the length of the femoral shaft in the GeomAge, MechAge, and MechGeomAge sets was adjusted according to Dragomir-Daescu et al. (2011). Nodes at the distal of the shaft were constrained to an external revolute joint in the coronal plane (Altai et al., 2019; Figure 2). Finally, a displacement rate of 100 mm/s was applied to the femoral head using a cement pad. The cement pad was free to move in the vertical direction, according to Dragomir-Daescu et al. (2011), while a fully constrained pad was added to the trochanteric area of the femur. Later, the force-displacement curves and the corresponding fracture regions were compared to the experimental results (Courtney et al., 1995;Dragomir-Daescu et al., 2011) in order to verify the results. The maximum force in  2000) and Mayhew et al. (2005), (Angle of rotation was assumed in such a way that the posterior of the femoral neck corresponds to 270°), (C) transverse plane of the shaft and neck (Plane's normal is in red), the femoral head was assumed to be isotropic.
Frontiers in Mechanical Engineering | www.frontiersin.org May 2021 | Volume 7 | Article 691171 the force-displacement curve, which will be denoted as fracture forces, was compared to demonstrate the effect of each of the changes separately and combined. For additional verification, the age-based fracture forces from similar experimental setups were pooled to compare with the results of the MechGeomAge sets. All simulations were done using LS-Dyna 971 R8.0.0 (LSTC 2016), and each sideways fall simulation took about 7 h to complete on four cores.

RESULTS
The geometrical changes due to aging in the femoral neck affect the fracture forces differently for men and women (Figure 3). Combination of subperiosteal expansion and endocortical resorption caused a slight increase in fracture forces in men, whereas it caused a decrease in women (Figure 3). The fracture forces increased by 25 N per decade, 0.5% of the baseline, for men, while it decreased by 116 N per decade, 3% of the baseline, for women. The changes in mechanical properties due to aging progressively reduced the fracture forces in both men and women ( Figure 4). However, the fracture forces reduced differently in men compared to women. While the fracture force decreased 354.5 N per decade, 7% of the baseline, for men, it decreased 225.4 N, 5.7% of the baseline, for Women. It is noteworthy that the geometry was kept intact for various ages in the MechAge set.
Similar to the MechAge, the fracture forces decreased with aging for both men and women in the MechGeomAge set ( Figure 5). In general, women had a lower fracture force comparing to men of a similar age, whereas both men and women experienced a similar trend in the age-dependent fracture forces. The fracture forces decreased by 373.3 N per decade, 7.2% of the baseline, and 368.1 N per decade, 10.7% of the baseline, for men and women, respectively.
Force-displacement curves from the simulations of the MechGeomAge set were compared to the experimental results (Courtney et al., 1995;Dragomir-Daescu et al., 2011). Several key characteristics of the force-displacement curves were similar for the experimental results and simulations ( Figure 6A). The general shape of the simulation curves was similar to the experiments even though there were apparent differences between each of the simulations and the experimental results ( Figure 6A). Apart from the similarities in the forcedisplacement curves, the resultant age-dependent fracture forces remained within the experimental fracture forces (Figures 6B,C). Furthermore, the fracture region was also similar between the experiment  and simulations. The maximum force in the simulation curves typically corresponded to a fracture in the superior and superoanterior sections of the femoral neck. The fracture pattern was similar to those reported in the experiment Figure 7). This fracture region varied slightly between the sexes and various ages ( Figure 6).

DISCUSSIONS
Aging causes two major changes in the femoral neck, geometrical and mechanical properties. Numerical simulations made it possible to split the effect of the geometrical and mechanical properties due to aging. In this study, we conducted a distinct series of simulations to investigate the potential effect of each of those changes separately and combined. In order to account for the asymmetry of the yield and ultimate strains as well as Young's and shear modulus, we implemented a transverse isotropic continuum damage model for the cortical bone. Overall, the fracture forces decreased with aging when accounting for both geometrical and mechanical properties changes in men and women. Likewise, the fracture forces reduced with aging in men and women when accounting for only the mechanical properties changes. However, the geometrical changes caused a diverging effect on the fracture forces in each of the sexes with a slight increase in men and a decrease in women with aging. The geometrical changes in the femoral neck consist of the changes in cortical thickness, cortical cross-sectional area, and moment of inertia. Ito et al. (2011) investigated the changes happening in a longitudinal two-year cohort study and found a significant reduction of cortical cross-sectional area and moment of inertia in women. Manske et al. (2006) showed that the failure load is associated with the cortical crosssectional area and moment of inertia. Although cortical cross-sectional area and moment of inertia were not directly changed in this study, these changes were indirectly altered. The subperiosteal expansion increases the moment of inertia, which acts in favor of increasing the bone's fracture force. Conversely, the endocortical resorption decreases the cortical thickness, which leads to a thinner loading area and, consequently, a lower fracture force (Duan et al., 2003;de Bakker et al., 2009). Thus, these two geometrical changes are essentially working against each other. The results of the GeomAge set indicated that the subperiosteal expansion is the dominant factor in men, and it assists the bone to maintain its strength during aging; whereas the rate of reduction of the cortical thickness dominates the changes in the neck region in women; and it causes a decrease in the fracture forces. Moreover, women experienced lower fracture forces compared to men at a similar age when accounting for agedependent geometrical changes. This difference between the fracture forces increased in the later years. A possible explanation for these differences might be that the cortical thickness is changing dissimilarly in the corresponding sectors among men and women (Mayhew et al., 2005;Johannesdottir et al., 2013). This inhomogeneous change leads to a thicker cortical thickness for men compared to women, which itself increases the load-bearing capacity of the bone. Additionally, the widening of the neck width amplifies the effect of greater cortical thickness in men.
Despite the geometrical changes, the variations in the mechanical properties of the cortical and trabecular bone due to aging affect men and women in an analogous way (Seeman 2001;Morgan and Bouxsein 2008;Keaveny et al., 2010). Decreases in the modulus and ultimate strain make the bone more brittle and reduce the toughness of compact bone (Burstein et al., 1976;Martin 1993;Zioupos and Currey 1998;Villette and Phillips 2018), which could explain the reduction of the fracture forces. Another important finding is that the size of the bone has a key role in the reduction rate of fracture forces due to aging. In other words, the fracture forces decreased 5.7% from the baseline for women as a result of the mechanical properties changes due to aging, compared to men, for which it decreased 7% of the baseline.
Accounting for both mechanical and geometrical changes due to aging, the fracture forces reduced for both men and women. Women experienced a lower fracture force (average of all age groups: 3359 N) and a higher reduction rate (10.7% of the baseline) compared to men (4909 N, 7.2%). This finding is consistent with previous research, which found a smaller fracture force for women compared to men in general (Kanis   , 1990;Rosen et al., 1999;Bauer et al., 2006;Dragomir-Daescu et al., 2011;Nasiri and Luo 2016;Nasiri Sarvi and Luo 2017;Dragomir-Daescu et al., 2018;Kleiven 2020). Rezaei and Dragomir-Daescu (2015), in an experimental study, suggested 640.3 N and 859.3 N per decade decrease in the fracture forces for men and women, respectively. Lang et al. (2012), in a 5-years longitudinal study, found 122 N and 164 N decrease per 5 years in fracture forces for men and women, using FE modeling. Similarly, Keaveny et al. (2010) predicted 547 N and 607 N decrease per decade in the fracture forces for men and women (10 and 13% decrease of the baseline). The absolute decrease in the fracture forces in the current study was 373.3 N for men and 368.1 N for women. Despite the differences in these absolute values, the percentage of change per decade is in agreement with those studies.
In the current study, THUMS v4.02 was used to acquire the geometry of the femur. In the baseline model, we aimed to modify the primary model to address the shortcomings of the current model. Afterward, we investigated the potential effects of the major aging changes of the geometry and mechanical properties on the fracture forces of the femur. Meanwhile, the findings of the current study are subject to several limitations. First, the model was simplified at the structural and material properties level, although several previous studies were focused on sophisticated CT image-based modeling of the femur. At the structural level, the accuracy of the cortical thickness estimation could be improved using a recently developed cortical bone mapping (CBM) technique (Treece and Gee, 2015;Treece and Gee, 2018;Schileo et al., 2020). Similarly, it can be used to follow the age-dependent changes to the cortical bone with a higher resolution. At the material properties level, the mechanical properties were applied homogenously to the bone even though the CT images have shown the inhomogeneity of those properties. A similar approach of CT-based models could be used to assign heterogeneous mechanical properties to the cortical and trabecular bones using bone mineral density or porosity (Bayraktar et al., 2004;Hansen et al., 2008;Ariza et al., 2015;Grassi et al., 2016;Enns-Bray et al., 2018). It has been shown that the ultimate strain can benefit from a similar density-based power law , although numerous CT-based models have implemented the uniform asymmetric strain limits and have shown accurate results when comparing to the experiments (Bayraktar et al., 2004;Grassi et al., 2016). Despite the benefits of a higher resolution femur model and the improved accuracy of these models, this approach could be computationally expensive for implementing in the whole-body models. The fracture criteria for the damage model of the current study were asymmetric, transverse isotropic, and uncoupled. A possible improvement in the damage model could be made using the continuum damage model proposed by Hambli (2013). The proposed damage model was asymmetric and coupled; however, it lacks the transverse isotropic/orthotropic material properties and fracture criteria. Another limitation of this work is that the interindividual differences such as neck-shaft angle and neck length, which could play a role in the femoral strength, were not studied. A recent study by Dragomir-Daescu et al. (2018) suggested that at least neck-shaft angle is not improving the prediction of the femoral strength as much as sex, age and BMD do. The simplifications of the current study, especially in the homogenous material properties of the bone, could adversely affect the slope of the force-displacement curves in Figure 6A, and the limited range of the fracture forces in Figures 6B,C. Also, the smaller difference between the slopes of the most compliant and stiffest models compared to the experiments can have been affected by using average geometry for each age group instead of subject-specific geometries. Part of the wide range of the fracture forces for each decade of age in Figure 6B for men and Figure 6C for women corresponds to the quality of the tested bone such as T-score or BMD. The age in the current study refers to a healthy adult in each decade of age, and a wider range could be achieved for each age group if the mechanical properties of different age group were assigned to all GeomAge models. Finally, having a larger number of experimental force-displacement curves would offer a better validation of the model. Despite these limitations, the current model uses minor modifications to the primary THUMS v4.02 model, which enables the direct application of it in a whole-body model. Moreover, the simplicity of the model does not interfere with the accuracy.