1 Introduction

Bone fractures are common childhood injuries which have been estimated to account for 25% of all paediatric injuries (Cooper et al. 2004). Recent evidence suggests that childhood musculoskeletal injuries may be related to higher risk of osteoarthritis later on in life (Antony et al. 2016). Children aged 3 years and younger have limited ability to protect themselves. This combined with limited communication skills, makes it difficult to diagnose inflicted injuries. In fact, the majority of fractures due to inflicted injury occur in children younger than 2 years old (Carty 1997; Loder et al. 2006). This is a big social problem with serious consequences for affected children and their families (Jayakumar et al. 2010). Unfortunately, there has been little quantitative data on bone development or strength in this age range. Previous studies were limited by sample size and available techniques to investigate the mechanical response of paediatric bone. Most of these studies focused on evaluating the strength of bone through mechanical testing (e.g. bending), which is destructive. There is also a complete lack of information on infant’s bone strength between new born and 1-year-old. This age group is particularly challenging to characterise due to the rapid growth phase (the mass of a typical neonate more than doubles in the first year of life), accompanied by marked changes in anatomy and function. Therefore, research in bone strength for very young children (0–3 years old) is highly desirable using the latest engineering techniques. This quantitative data can be used to enhance our understanding of the injury tolerance of young bones under different loading conditions, which are frequently presented in both accidental and inflicted injuries. This can then be compared with the force resultant from the physical events described by the parents/carer, combined with a dynamic model in the future. Such data will also aid child dummy designs and the design of car restraint systems and create more accurate computer models of car crash tests in the future.

Fracture type is dictated by the direction and magnitude of the applied load. Long bones may fail in tension, compression, bending or torsion, and under a combination of all or some of these loading mechanisms. In case of an injury, the bone will experience a combination of loads due to structural asymmetry as well as uneven load transfer depending on the injury scenario. The latter is more complex and requires the construction of a full-body dynamic model to quantify the amount of force acting on the femur for each injury scenario. This could be investigated in the future using open source child models such as PIPER (http://www.piper-project.eu/), but is beyond the scope of this study. The former can be addressed to some extent by separating these loading conditions and testing each of them independently in order to quantify a range of failure load, which is the approach taken in this study. For example, a transverse fracture (characterised by a fracture line perpendicular to the long axis of the bone) may occur as a result of tensile loading and bending, whereas a torque applied to the femur may lead to spiral fracture (Gitajn and Rodriguez 2011; Pierce et al. 2004; Turner and Burr 1993). Therefore, by investigating the bending and torsional strength of paediatric femur, the results of this study will provide an indicative failure range under these loads, either separately or combined.

Previous studies have frequently used impact or quasi-static loading on the bone in order to mimic these loading conditions in a well-controlled experimental set-up. To date, only a few studies have investigated fracture tolerance of immature human bones using an experimental approach (Forman et al. 2012; Miltner and Kallieris 1989; Ouyang et al. 2003). Miltner and Kallieris (1989) conducted three-point bending tests, using quasi-static loading at a rate of 50 mm/min on 28 cadaveric lower limbs (bones with surrounding soft tissues). Samples were taken from children aged 1 day to 6 years old. The reported fracture moment ranged from 7.05 to 109.5 Nm. In 2003, Ouyang et al. performed three-point bending tests under quasi-static loading on eleven pairs of cadaveric humeri, radii, ulnae, femora, tibiae, fibulae of children aged 2–12 years old. Bones were freshly isolated by removing the surrounding soft tissues. The maximum bending force reported for the femora, tibiae and fibulae was in the range of 15.7–72.6, 11.3–51.3 and 2.4–20.8 Nm, respectively. The most recent experimental study on paediatric bone strength was reported by Forman et al. (2012), in which three-point bending tests were conducted on 23 cadaveric femora aged between 1 and 57 years, using a loading rate of 1.5m/s. They reported a bending fracture moment ranging from 61.4 (children) to 675 Nm (adult).

With recent developments in advanced medical imaging and in silico technologies, non-invasive techniques such as computed tomography (CT)-based finite element (FE) analysis have been extensively applied to study bone strength in adults, particularly in the application of fall prediction (Grassi et al. 2012; Lotz et al. 1991; Lotz and Hayes 1990). The method has been fully validated against experimental data collected from adult bones (Mccalden et al. 1997). These studies were built on a positive correlation between bone density and modulus of elasticity. This relationship has recently been confirmed in children in a study conducted by Öhman et al. (2011), these authors suggested, through an experimental study on human bone tissues, that the correlation between ash density, both strength and stiffness found in adults can be extended and applied to children’s bones. Based on this evidence, a previous study conducted by Li et al. (2015) reported the development of a framework to model paediatric long bones and showed initial bending analysis of 15 femora (0–3years old) using a CT-based FE approach (Taddei et al. 2004; Schileo et al. 2008).

During the same period, there were a few other FE studies on paediatric bones, with one on infants and the others on older children (6 years and above). Yadav et al. (2017) investigated the effect of muscle groups’ activation on the growth of femur using MRI-based FE models of femora for three children (aged 6, 7 and 11 years). Both Meng et al. (2017) and Angadi et al. 2015 investigated the accidental injury of older children (e.g. pedestrian accidents). The former used scaled models, generated from adult data, while the latter used paediatric bone models made of simplified composite materials. Scaled-down geometries from adults are not representative of children as they have very different anatomy. Furthermore, it is evident from the literature that there is an absence of personalised anatomical or material data on infants and very young children during the rapid growth phase, when they are most vulnerable to inflicted injuries.

Therefore, the primary aim of this study is to quantify the bone strengths (under bending and torsion) for 30 children aged 0–3 years old using a CT-based FE approach. The personalised mechanical properties will be estimated from CT attenuation. A preliminary relationship between bone strength and age/mass will be determined from the simulation results. The information reported here will inform future FE studies on infant and paediatric bones. These models may also be used to create surrogate models of very young children in the future, and combined with dynamic models to simulate different injury scenarios.

Table 1 Demographics for the current cohort, consisting of 30 cases

2 Material and methods

2.1 Finite element model generation

Post-mortem CT scans of 30 children performed at the Sheffield Children’s Hospital were used for this study. This data set has been expanded since previously reported by Li et al. (2015). Table 1 shows the geometry of the femur (length and cross-sectional area), and the material properties (modulus of elasticity) in the current cohort.

The right femur of each child was segmented using ITK-Snap (ITK, http://www.itksnap.org). The segmented bone surface was meshed with 10-node tetrahedral elements in ICEM CFD 15.0 (Ansys INC., PA, USA). Material properties were mapped from the CT scans using a well-established material-mapping procedure (Bonemat v3, Rizzoli Institute) (Schileo et al. 2008). The peak elastic modulus value in the cohort ranged between 14 and 22.5GPa. Detailed methodology is described in Li et al. (2015).

Fig. 1
figure 1

Reference system used to align the finite element model for bending and torsion (anterior view); a The original reference system reported in Li et al. (2015); the x-axis passes through the long axis of the distal half of the femoral shaft but misaligns in the proximal half of the shaft. b The improved reference system where the x-axis is now defined as a line passing through the centroids of two cross-sectional areas at 25% and 75% of total femoral length. Z-axis is perpendicular to the x-y plane and points anteriorly (out of the page) following the right-hand rule

2.2 Reference system

Due to the asymmetry of the bone, a perfect alignment of the coordinate system with the long axis of the shaft would be difficult to achieve. The previous reference system reported in Li et al. (2015) ensured a good alignment in the distal half of the femoral shaft, with less accuracy towards the proximal end (Fig. 1a). Here an improved coordinate system was proposed, where two cross sections of the femoral shaft were identified at 25 and 75% of the total length (distance between the proximal and distal ossification centres). The centroids of these two cross sections were estimated and the x-axis (long axis) was defined to be a line passing through both centroids. This will ensure better alignment of the x-axis with the major axis of the shaft in an attempt to minimise the bending effect due to geometrical asymmetry. Positive y-axis points medially, from the mid-shaft cross section to the proximal ossification centre. Positive z-axis points to the anterior, which is perpendicular to the xy plane (Fig. 1b). The proposed improved alignment system is similar to an alignment approach previously optimised in (Cheong and Bull 2015) in order to align long bones in an experimental testing of pure bending and pure torsion. This is different from the conventional reference system used for adults (Hilal and Leardini 1999), due to the lack of features in the epiphysis of infant’s femurs. Consequently, a different set of landmarks were required in order to set up a reliable reference system unique to children’s anatomy.

2.3 Boundary conditions

In order to investigate the effect of bending and torsion, two different boundary conditions were defined. Four-point bending was simulated by applying two equal forces to the femoral shaft in the y direction where the span of the loading is equal to half span of the supports (see Fig. 2a). Various orientations around the longitudinal axis of the shaft were then analysed with a \(10^{{\circ }}\) increment (35 increments in total). Only the femoral shaft (mineralised portion of the bone) was used for this simulation, which represented approximately 50% of the total length of the femur (Li et al. 2015). This was because the proximal and the distal ends of the paediatric femur were largely composed of materials in the transition to fully mineralised bone. The contribution to bending strength from non-mineralised bone would be much lower compared with the mineralised region. Furthermore, under the current boundary condition, this part of the bone would appear redundant in the simulation and was removed from the analysis.

Fig. 2
figure 2

Schematic of the bending and torsion loading conditions; a Four-point bending: two equal forces applied to the femoral shaft, adapted from Li et al. (2015); and b Torsion: moment around the x axis was applied at the distal end while the rotational movement of the proximal end was fixed in this direction. Regions of interest for both bending and torsion are highlighted in the red boxes

Fig. 3
figure 3

First principal strain distributions under bending and torsional loads are shown in the lateral and medial views. The illustrated femur model was obtained from a 1-year-old. The direction of loading was applied on the medial side of the femur in the bending simulation. External rotation was applied for the torsional simulation

For torsion, the full length of the femur was used in order to set up an adequate boundary condition. Two pilot nodes were added to the finite element model to define the axis of rotation along the femoral shaft. The pilot nodes were defined at the proximal and distal ends of the femur along the x-axis. The nodes were created at a distance equal to half of the total length between the centroids (at 25 and 75% of total length) (see Fig. 2b). A multi-point constraints (MPC) method was used to connect the distal and proximal nodes of the femur with the pilot nodes. The linear constrained equation related all the degrees of freedom (translation and rotation in x, y, and z directions) of nodes at both ends to the pilot nodes. The distal pilot node was fixed along the longitudinal translation (x-axis) and all other degrees of freedom were free. The proximal pilot node was constrained in y and z directions, and was not allowed to rotate in the direction of the applied moment. A torsional moment was applied through the distal pilot node with two different rotational directions representing internal and external torsion of the child’s leg.

2.4 Post-processing

For bending, the force to failure was calculated following maximum principal strain criterion for each orientation. As suggested by previous studies (Schileo et al. 2007), failure of the bone occurs when the maximum strain reaches the threshold of the elastic strain limit. For human bone these limits have been reported to be 0.73% in tension and 1.04% in compression (Bayraktar et al. 2004). For comparison purposes, the equivalent failure moment for bending was also calculated.

For torsion, the moment to fail was estimated for both internal and external torsion. The locations of maximum tensile and compressive principal strains under the two different rotational directions were recorded for each case. Regions of interest are highlighted in Fig. 2.

3 Result

For both bending and torsional loading, the maximum first principal strain was observed at the mid-shaft, where the cross-sectional area of the shaft is the narrowest. The distribution of the first principal strain is illustrated in Fig. 3 using the 1-year-old femur model. For all 30 cases, the moment to failure under bending ranged from 0.85 to 27.9 Nm, with an equivalent load to failure of 97–1022 N. The moment to failure for external and internal rotation was very similar, with 1–31.4 and 1–30.7 Nm, respectively. For the majority of the cases, the predicted load and moment to failure increased consistently with age, apart from one outlier, illustrated with a cross in Fig. 4. This outlier was 40 weeks old with noticeably lower body mass (\(< 9\hbox {th}\) percentile) and height (\(<2\hbox {nd}\) percentile), and a corresponding low peak modulus (as shown in Table 1). The data suggest that this individual may experience some growth delay. Plotting the above results against the mass removes this outlier.

In general, the load and moment to failure increased with age and mass. However, because of the limited number of cases and the non-uniform distribution of the age within the cohort (more cases under 1 year old), it was difficult to conduct a robust regression analysis. Therefore, the resulting data were fitted with two possible regressions: one linear and one quadratic, as shown in Fig. 4. The same procedure was then repeated using only data from samples between 0 and 6 month old, in order to show the trend in the infant age range, as shown in Fig. 5. The moment to failure under bending predicted by this study agreed reasonably well with previous experimental data reported by Ouyang et al. (2003) and Forman et al. (2012), given the small number of cases reported between 1 and 3 years old (see Fig. 6).

Fig. 4
figure 4

Predicted moment to failure (Nm) plotted against the age and mass for bending and torsion; two possible regressions were indicated: linear and quadratic. The outlier is indicated with a cross. The torsion results shown are for external rotation; however, a similar trend was observed for internal rotation. The p value of all regressions are \(<\,0.001\)

4 Discussion

CT-based finite element models of paediatric femora were developed and simulated under bending and torsional loading in order to provide personalised information on bone strength for infants and very young children during the rapid growth phase. The results showed good agreement between the current predicted values and the previously measured experimental values on children’s bones reported in the literature (Forman et al. 2012; Miltner and Kallieris 1989; Ouyang et al. 2003). The findings add invaluable quantitative data of femora strength to the scarce literature to date.

Both cross-sectional area and length of the femur increased with age/mass (as shown in Table 1). This increase in femur size is likely to increase the bending and torsional strengths, as shown by the simulation results. However, it was observed that developmental level was heterogeneous in the current cohort. Consequently, body mass appeared to be a better indicator for bone strength, compared with age. The results were also suggestive of a faster increase in bone strength beyond 1 year of age (this needs to be confirmed with further cases). It is well known that mechanical factors are strongly associated with bone growth. The increase in magnitude and frequency of loading, as well as the increase in mobility in the early years, would strongly encourage bone growth in order to adapt to the evolving load. This is reflected by the increase in femur length and cross-sectional area (as shown in Table 1), accompanied by continuous mineralisation of bone tissues.

Fig. 5
figure 5

Moment to failure (Nm) plotted against the age and mass under bending and torsion for age ranged between 0 and 6 months old to illustrate the trend in the infant age range. The p value of all regressions are \(< 0.001\)

The predicted moment to failure of the current cohort, which represents bone strength, increased steadily with age and mass. As shown in Figs. 4 and 5, the predicted trend was in agreement with previous experimental studies, where a positive relation between fracture moment and age was also reported (Forman et al. 2012; Ouyang et al. 2003). The range of fracture moment predicted by the current FE model under four-point bending was 0.85–27.9 Nm. Although this range is comparable to that reported in Ouyang et al. (2003) and Forman et al. (2012) (see Fig. 6), it should be noted that the bending moments reported in these two studies were measured under three-point bending tests, which were likely to predict higher strength than four-point bending. Three cases from Ouyang et al. (2003) were within the age range of the current study, at 2, 2.5 and 3 years, respectively. The reported fracture moments were 29.6, 24.3 and 39.6 Nm, respectively. Three cases reported in Forman et al. (2012) fell within our age range, one at 1.33 years old and another two at 2 years old. The fracture moments of these femora were 61.4, 61.7 and 65.5 Nm, respectively. The reason for the comparably higher fracture moment predicted by Forman’s study could be due to the much faster loading rate used (1.5 m/s compared with 0.008 m/s in Ouyang’s study). Furthermore, Ouyang’s study was based on a Chinese population, whereas Forman’s study was based on a Spanish population. Consequently, ethnicity may play a role in the measured differences in bone strength. It should be noted that there was no comparable experimental data for children below 1 year old, where the majority of the current cohort sits. To the authors’ knowledge, there is only one study in the literature that reported three-point bending tests on 28 children aged between 1 day and 6 years old (Miltner and Kallieris 1989). Quasi-static loadings were used for 18 specimens while dynamic loadings were used for the remaining 10 tests. This study reported a fracture moment ranged from 7.05 Nm (6 days) to 109.5 Nm (6 years). However, the entire lower limb was used in this study including all surrounding soft tissues. Previous literature suggests that soft tissues (including skin, fat and muscles) surrounding the bone may absorb some energy during loading. Therefore, the results of Miltner’s study are not readily comparable with those reported here (Kerrigan et al. 2003)

Fig. 6
figure 6

Moment to failure for bending as a function of age predicted by the finite element model of the current study, plotted together with the measured moment to failure by previous experimental studies

We could not find any study in the literature that performed a torsional loading test on paediatric bones in order to compare with the current study. However, the predicted failure moment under torsion in the current study was consistent compared to the bending results, showing both a similar trend with age/mass and the range of predicted moments to fail. The predicted maximum principal strain under torsion was located at the narrowest cross-sectional area of the diaphysis, namely the mid-shaft. This is in agreement with previous experimental studies of long bones under pure torsional loads (Kress et al. 1995; Theobald et al. 2012).

The current study has a number of limitations. One of these is the non-uniform distribution of the cohort’s age. There were only four cases between 6 months and 3 years. This has prevented a robust regression analysis to be conducted on the relationship between bone strength and age/mass. Two possible relationships were suggested here: linear and quadratic. However, it is difficult to determine a starting point for failure moment at birth (or at zero age). One would imagine that the bone strength should be nonzero at birth. This renders some of the linear regression results invalid, as they produced negative bone strengths at this age. This needs to be considered when interpreting the regression results in the current study. In the future, more cases will be added to the current cohort to further enhance the predictive accuracy of the regression analysis.

Another limitation is related to the selected failure criterion. It was assumed that the cortical bone in children fails under peak strain in a similar manner as adult, and therefore a similar failure strain criterion (used in adults) could be applied to children. Although the peak elastic modulus value in infants is comparable to adults (Li et al. 2015), It should be noted that the assumption behind the failure criterion used (bone is a fragile material) is somehow conservative with respect to the amount of energy that bone can absorb before macroscopically failing. If in reality children’s bones fail in a more ductile manner, the total load would be higher, but if failure was assumed when the proportionality limit is exceeded, the results would be quite similar to those reported here.

Due to scarcity of data, it was assumed that the high correlation between the ash density and both strength and stiffness found in Öhman et al. (2011) (4–15 years old) can be extrapolated to represent the age range in the current study (0–3 years old). The predicted failure strength was then compared against previous experimental data (Forman et al. 2012; Miltner and Kallieris 1989; Ouyang et al. 2003), which were confirmed to be within the same range. Future work will need to be conducted in order to confirmation the validity of this assumption and validate the FE model. One potential approach could be through the testing of young animal bone, as reported in Cheong et al. (2017).

The current study suggests that the finite element approach, which has been widely and successfully used in adults, can be adapted and applied to children. Since getting paediatric bone samples is very difficult, using non-invasive techniques such as CT-based finite element analysis provides a valuable alternative to the investigation of paediatric bone biomechanics. In future, this technique will allow us to create surrogate models for infants and very young children, as well as obtaining more quantitative information about bone growth and strength to drastically enhance the little information currently available in the literature.