Martin Revel

Marc Gardegaront

François Bermond

David Mitton

Hélène Follet
email: helene.follet@inserm.fr

A credible homogenized finite element

Keywords: Forward Fall, Finite Element Analysis, Fracture Prediction, Radius, Bone strength Highlights

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Graphical Abstract

1 Introduction

Fragility fractures represent a worldwide health problem in the elderly (65 years and older) and remain underdiagnosed and untreated [START_REF] Khosla | A Crisis in the Treatment of Osteoporosis[END_REF][START_REF] Roux | Osteoporosis in 2017: Addressing the crisis in the treatment of osteoporosis[END_REF]. Occurring after a fall from standing height or lower (considered as low trauma fractures), radius fractures are the most common upper extremity fractures [START_REF] Vogt | Distal radius fractures in older women: a 10-year follow-up study of descriptive characteristics and risk factors. The study of osteoporotic fractures[END_REF]. These fractures are associated with age-related decline in bone quality, and because they occur earlier in life than other osteoporotic fractures, they can be interpreted as a warning signal for later and more deleterious fractures [START_REF] Melton | Assessing forearm fracture risk in postmenopausal women[END_REF].

Subject-specific finite elements models have been proposed to evaluate bone strength and strain, and fracture risk. Currently, the most widely used tool in clinical research is High Resolution peripheral Quantitative Computed Tomography (HR-pQCT). Based on a restricted area of approximately 9 mm length of the radius (located 9.5 mm to 18.5 mm proximal to the endplate of the radius), micro-finite element (µFE) analysis [START_REF] Boutroy | Finite element analysis based on in vivo HR-pQCT images of the distal radius is associated with wrist fracture in postmenopausal women[END_REF][START_REF] Whittier | Harmonizing finite element modelling for non-invasive strength estimation by high-resolution peripheral quantitative computed tomography[END_REF][START_REF] Whittier | Guidelines for the assessment of bone density and microarchitecture in vivo using high-resolution peripheral quantitative computed tomography[END_REF][START_REF] Vilayphiou | Finite element analysis performed on radius and tibia HR-pQCT images and fragility fractures at all sites in postmenopausal women[END_REF][START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF][START_REF] Arias-Moreno | Validation of distal radius failure load predictions by homogenized-and micro-finite element analyses based on second-generation high-resolution peripheral quantitative CT images[END_REF] can be performed to assess quasi-static axial bone strength. These µFE models showed good prediction (R² between 0.73 and 0.92) when comparing numerical outcomes with quasi-static experimental data [START_REF] Van Rietbergen | A survey of micro-finite element analysis for clinical assessment of bone strength: the first decade[END_REF]. However, despite this good level of prediction, retrospective [START_REF] Melton | Assessing forearm fracture risk in postmenopausal women[END_REF][START_REF] Boutroy | Finite element analysis based on in vivo HR-pQCT images of the distal radius is associated with wrist fracture in postmenopausal women[END_REF][START_REF] Vilayphiou | Finite element analysis performed on radius and tibia HR-pQCT images and fragility fractures at all sites in postmenopausal women[END_REF][START_REF] Chevalley | Fracture history of healthy premenopausal women is associated with a reduction of cortical microstructural components at the distal radius[END_REF][START_REF] Nishiyama | Women with previous fragility fractures can be classified based on bone microarchitecture and finite element analysis measured with HR-pQCT[END_REF][START_REF] Vilayphiou | Finite element analysis performed on radius and tibia HR-pQCT images and fragility fractures at all sites in men[END_REF] and prospective [START_REF] Sornay-Rendu | Bone Microarchitecture Assessed by HR-pQCT as Predictor of Fracture Risk in Postmenopausal Women: The OFELY Study[END_REF][START_REF] Ohlsson | Cortical Bone Area Predicts Incident Fractures Independently of Areal Bone Mineral Density in Older Men[END_REF][START_REF] Biver | Evaluation of Radius Microstructure and Areal Bone Mineral Density Improves Fracture Prediction in Postmenopausal Women[END_REF] cohort studies have not been able to establish better prediction of fracture risk with µFE analysis compared to the standard reference for the detection of osteoporosis fracture, which is the measurement of areal Bone Mineral Density (aBMD) by Dual X-ray Absorptiometry (DXA). It is only with the BoMIC consortium study, including participants of eight cohorts, that µFE outcomes associated with femoral neck aBMD improved the ability to predict non-traumatic or traumatic incident fracture by enhancing the area under the receiving operating characteristic (ROC) curve (p < 0.05) [START_REF] Samelson | Cortical and trabecular bone microarchitecture as an independent predictor of incident fracture risk in older women and men in the Bone Microarchitecture International Consortium (BoMIC): a prospective study[END_REF]. Nevertheless, this was not the case when considering only major osteoporotic fractures (low-trauma fractures) or when compared µFE outcomes with ultra-distal radius aBMD [START_REF] Samelson | Cortical and trabecular bone microarchitecture as an independent predictor of incident fracture risk in older women and men in the Bone Microarchitecture International Consortium (BoMIC): a prospective study[END_REF].

Radius bone strength through axial quasi-static load is associated with incident fractures, but despite perfectly predicted bone strength, fractures cannot be accurately predicted [START_REF] Van Rietbergen | A survey of micro-finite element analysis for clinical assessment of bone strength: the first decade[END_REF], since fractures in the distal radius mainly result from a fall from a standing height in a forward direction [START_REF] Oskam | Fracture of the distal forearm: epidemiological developments in the period 1971-1995[END_REF][START_REF] O'neill | Age and sex influences on fall characteristics[END_REF].

Our recent study examined the influence of loading conditions (angle of compression or torsion) in µFE models of 9.02 mm of the distal radius [START_REF] Revel | Influence of Loading Conditions in Finite Element Analysis Assessed by HR-pQCT on Ex Vivo Fracture Prediction[END_REF]. The results suggested that this approach does not provide improvement in classifying ex vivo fractured and non-fractured radii within an elderly population. Using soft layers at the cut faces of the HR-pQCT scanned section, Christen et al. [START_REF] Christen | Subjectspecific bone loading estimation in the human distal radius[END_REF] showed that axial loading on the wrist induced axial and shear forces as well as torsion and bending moment on the scanner's region of interest. These authors also showed that physiological bone loading conditions on the cross-section were patient-specific and need a full distal radius scanned region for the determination of realistic load transfer to the scanned region, which is not currently feasible in vivo with HR-pQCT. For this reason, we suggest that enlarging the scanning region may facilitate observation of loading sensitivity and reproduction of the in vivo mechanism occurring during a fall that led to fracture.

Other types of radius FE models have been proposed in a larger region of interest with a homogenized finite element analysis [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF][START_REF] Arias-Moreno | Validation of distal radius failure load predictions by homogenized-and micro-finite element analyses based on second-generation high-resolution peripheral quantitative CT images[END_REF][START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF][START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF]. The accuracy of these homogenized FE (hFE) models is heavily dependent on material properties, boundary conditions, and chosen failure criterion that define fracture threshold [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF]. However, their use have faded since the introduction of the HR-pQCT µFE during the decade of 2000-2010. Indeed, although hFE models have predicted bone stiffness, strain, or fracture strength obtained from static experimental tests, they have never been used to predict incident fracture occurring after an impact on the radius resulting from a forward fall on a given population.

Strain rate, in quasi-static constant displacement rate or impact loading, are important considerations for the mechanical behavior of the bone. Compact bovine bone femur was three to four times more brittle under a dynamic load than under a quasi-static load [START_REF] Pithioux | Comparison of compact bone failure under two different loading rates: experimental and modelling approaches[END_REF].

Moreover, impact velocity has an influence on trabecular bone strain rate which affects the compressive ultimate stress [START_REF] Enns-Bray | Strain rate dependency of bovine trabecular bone under impact loading at sideways fall velocity[END_REF]. On human femoral bone, relationships can be found between impact and quasi-static failure loads. An impact scenario can be extrapolated by quasi-static tests with errors of less than 12% [START_REF] Jazinizadeh | Comparing the fracture limits of the proximal femur under impact and quasi-static conditions in simulation of a sideways fall[END_REF]. We thus hypothesized that a quasi-static homogenized subject-specific finite element model of the distal radius could, with a different scanner modality than usually used on this site, predict bone fracture following a forward fall.

To confirm this hypothesis, this study aims to propose a credible finite element model to discriminate fractured and non-fractured radius groups obtained experimentally after a forward fall [START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF]. The credibility of the model was assessed regarding the surface strains field in specific areas and for failure load prediction. The model sensitivity was assessed in terms of the elasticity-density equation and failure criterion.

Materials and Methods

Specimens and Image Acquisition

Data were acquired from 30 fresh-frozen cadaveric left radii from elderly donors (50-96 y.o.; 15 males, 15 females) provided by the University Department of Anatomy Rockefeller (Lyon, France) and were considered for this study (French Ministry of Education and Research, authorization n°DC-2015-2357). The most distal portion of the bones cleaned of soft tissues were imaged with a clinical CT scanner (New-Tom 5G, QR, Verona, Italy). The New-Tom 5G device was chosen for its good spatial resolution in the spectrum of a dental cone beam CT with a nominal detector resolution of 127 µm x 127 µm. Images were acquired on a segment of 8 cm length in the longitudinal direction beginning at the distal styloid of the radius using the standard protocol proposed by the manufacturer, with a tube voltage of 110 kVp and adaptive mAs. Acquisition time was 24 s with an effective dose of approximately 10 μSv (mean radiation dose (CTDI) = 12.6 mGy). Images were reconstructed with a field of view of 8.16 cm and a matrix of 544×544 to give a reconstructed isotropic voxel size of 150 μm.

Scans of the thirty specimens were made along with a forearm calibration phantom (QRM, Moehrendorf, Germany), with calcium hydroxyapatite equivalent concentrations of 0 and 200 mg/cm 3 . These scans were used to establish the following linear relationship (i) between CT Hounsfield units (HU) and calcium hydroxyapatite equivalent density 𝜌 𝐻𝐴 in g/cm 3 :

𝜌 𝐻𝐴 = 0.00038 × 𝐻𝑈 -0.1262 (i)
Each radius was wrapped in saline-moistened gauze and frozen at -20°𝐶 before the experiments.

Ex-vivo Experiment

The experimental protocol used here was previously described in detail by Zapata et al. [START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF].

The day before the experiments took place, bones were thawed for 16 h at 4°C and then 6 h at room temperature. Using a positioning laser, radii were potted with an 8 cm free length in a polyurethane resin in a steel cylinder with an angle of 75±2° between the anterior face of the radius and the impactor. This position reproduces the most common orientation in a forward fall [START_REF] Chiu | Prediction of upper extremity impact forces during falls on the outstretched hand[END_REF]. A rigid polyurethane mold was made to reproduce a simplification of the radio-carpal joint, which is involved in the fracture mechanism [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF], in order to evenly distribute the load on the articular surface of the radius throughout the impact.

In order to limit the loading of the radius so as to avoid having bone fracture in all cases, the potted bones were placed in a cylinder bar on a rail system, which was free to slide along the loading axis. The weight of the bar (12.5 kg) was considered as an arbitrary weight to represent the mass involved during a fall. This mass was identical for all the tests. A six-axis sensor (105515TF, Humanetics, Germany) was tightened to the impactor located in front of the rail system. The radii were then loaded through the resin mold at 2 m/s using a hydraulic high-speed testing machine (LF technologies, France). The initial distance between the impactor and the mold was 50 mm which allowed the acceleration of the impactor and stabilization of its speed to reach 2 m/s before impact. The displacement of the stroke was set to 10 mm before stops.

Two laser scan acquisitions were made for each specimen. The first laser scan was made on the articular mold at the surface in contact with the articular surface of the radius. This procedure was performed before fixation of the articular mold on the radius using a silicon rubber (Figure 1-A). The second laser scan corresponds to the whole radius after positioning the pot into the cylinder, just before the impact (Figure 1-B). These two laser acquisitions were useful for the finite element model to retrieve specific loading surface on the radius and specific orientation at the beginning of the experiments.

Observations were recorded using four high-speed cameras (FASTCAM SA3, Photron, Japan) at 2 000 frames per second using a shutter speed of 50 µs with a resolution of 1024 x 1024 pixels (∼24 pixel/mm). The cameras were placed in pairs with an angle of 30° (Figure 1-C) and with a distance to specimen of approximately 540 mm. The front cameras, facing the ulnar surface of the radius, recorded the test using 105 mm F2.8 DG Macro sigma lenses. The upper cameras, facing the anterior surface of the radius, recorded the test using 50 mm Z1.4/50 mm ZF planar Zeiss lens. Three spotlights (400D, Dedolight, Germany) provided sufficient lighting to avoid any shadow areas on the images and to ensure good image quality.

A speckle was previously painted on the bones to compute strain fields on a specified nonfractured area of the anterior face of the radius (Figure 1-D) by analyzing the video recordings using VIC3D stereo-correlation software (Correlated Solution, South Carolina, USA). A calibration plate was used before each test to re-calibrate the systems for each specimen. The resulting average projection error was ∼0.02 pixels. DIC settings were chosen according to manufacturer recommended guidelines: facet size was chosen based on an average valued of ∼0.01 pixels of the correlation tracking function (sigma estimate) and resulted to ~41 pixels.

Step size (referred to overlap) was roughly a quarter of the subset and resulted to ~10 pixels. The noise through time between consecutive images was measured with respect to the reference image. Thus, a time smothering filter was applied (More detailed information about DIC-settings are provided in Supplementary Table 1).

All samples were impacted under the same loading conditions leading to two groups: fractured and non-fractured radii.

FE Model

Models of the radii based on segmented CT data from the open-source platform 3D Slicer [START_REF] Fedorov | 3D Slicer as an Image Computing Platform for the Quantitative Imaging Network[END_REF] were imported into Ansys Software (v. 2019 R1 Ansys® Inc., USA) to be meshed (10-nodes, tetrahedra) and to perform the numerical analysis. These continuum meshes were created respecting the quality recommendations on the distortion of the elements by inspecting three parameters: aspect ratio, angle idealization, and element Jacobians [START_REF] Burkhart | Finite element modeling mesh quality, energy balance and validation methods: a review with recommendations associated with the modeling of bone tissue[END_REF]. In order to determine the optimal number of elements for the resolution, a mesh convergence study was carried out on a sub-sample of 6 specimens (3 men and 3 women) ranging from 50 000 to 600 000 elements by increments of 50 000. Two parameters were evaluated: mean equivalent von Mises strain in a region of interest of the anterior face of the radius and the maximum von Mises stress in the whole model excluding edge effects (excluding the 5% of the highest stress values). The convergence of the model was considered when the output parameters showed less than a 2% difference between three successive iterations and finally resulting in 200 000 elements, leading to a mean-size element of 0.9 mm edge-length (Supplementary Material).

In order to evaluate the influence of elastic property on the FE models, three previously established relationships were investigated (Eqs. ( 1)-( 3)) to assign inhomogeneous elastic material properties to the FE models (Bonemat v3, [START_REF] Taddei | The material mapping strategy influences the accuracy of CT-based finite element models of bones: an evaluation against experimental measurements[END_REF]). Eq. ( 1) was developed by Keller et al. on lumbar vertebrae and femoral metaphyses and diaphyses [START_REF] Keller | Predicting the compressive mechanical behavior of bone[END_REF], Eq. ( 2) was developed by Morgan et al. on femoral neck [START_REF] Morgan | Trabecular bone modulus-density relationships depend on anatomic site[END_REF], and Eq. ( 3) was obtained by Duchemin et al. on human cortical femur [START_REF] Duchemin | Prediction of mechanical properties of cortical bone by quantitative computed tomography[END_REF].

𝐸 = 10 500𝜌 𝑎𝑠ℎ 2.57 [34], (1)
𝐸 = 6 850𝜌 𝑎𝑝𝑝 1.49 [35], (2)
𝐸 = 10 095𝜌 𝐻𝐴 [36], (3)
where Young's modulus 𝐸 is expressed in MPa, and the densities in g/cm 3 , 𝜌 𝑎𝑠ℎ (ash density), 𝜌 𝑎𝑝𝑝 (apparent density), and 𝜌 𝐻𝐴 (hydroxyapatite equivalent density, obtained with Eq. (i)).

Calcium hydroxyapatite equivalent density obtained from the calibration was converted to ash density using:

𝜌 𝑎𝑠ℎ = 0.0698 + 0.839 𝜌 𝐻𝐴 [START_REF] Kaneko | Mechanical properties, density and quantitative CT scan data of trabecular bone with and without metastases[END_REF].

Ash density was normalized to apparent density with a ratio of 0.6 between 𝜌 𝑎𝑠ℎ and 𝜌 𝑎𝑝𝑝 , within the range 0.55-0.63 identified in the literature [START_REF] Keyak | Correlations between orthogonal mechanical properties and density of trabecular bone: Use of different densitometric measures[END_REF][START_REF] Schileo | Subject-specific finite element models implementing a maximum principal strain criterion are able to estimate failure risk and fracture location on human femurs tested in vitro[END_REF]. This resulted in around 220 materials ranging from 0.05 to 12.1 GPa for Eq. ( 1), from 0.01 to 15.7 GPa for Eq. ( 2), and from 0.03 to 12.0 GPa for Eq. (3) (Figure 2). Each element was assigned a Poisson's ratio of 0.3.

The models were specifically oriented in the same configuration as the experiments (with an uncertainty of ±2°) by matching surfaces between the laser scan acquisition and the FE model. The specific maximal peak load observed experimentally was applied to the articular surface along the z axis and the proximal end of the radius at the location of potting was fully constrained (Figure 3). Analysis were performed into Ansys Software (v. 2019 R1 Ansys® Inc., USA).

Three numerical failure criteria were evaluated in this study in order to assess their sensitivity to discriminate fractured bones from non-fractured bones. These criteria have shown good prediction of quasi-static radius [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF] and femoral [START_REF] Duchemin | Prediction of mechanical properties of cortical bone by quantitative computed tomography[END_REF] failure load. An element failure occurred when a criterion is greater than or equal to 1 (Table 1 Failure loads were then back calculated, assuming that fracture occurred when a cluster of contiguous yielded elements exceeded a predefined volume. Different failed volume has been proposed in the literature: approximatively 150 𝑚𝑚 3 for micro finite element models [START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF] to 350 𝑚𝑚 3 for homogenized FE models [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF][START_REF] Troy | Off-axis loads cause failure of the distal radius at lower magnitudes than axial loads: a finite element analysis[END_REF] of the distal radius. A variation of defective volumes ranging from 150 to 450 𝑚𝑚 3 in increments of 100 𝑚𝑚 3 was thus examined. The numerical failure forces obtained for each case was then integrated into the sensitivity study for the detection of fractures.

Data Analysis

Surface Strain Comparison

Experimental surface strain field was evaluated at the experimental measured peak force for non-fractured bones. The areas from the digital image correlation (DIC) were matched with the surface of the FE models by a custom iterative closest point approach developed in our laboratory in order to compare the bone numerical response with the experimental data analyzed by stereo-correlation. Principal surface strains measured at peak force were retrieved. Points where the confidence interval σ (a measure for how accurately the displacement of each point is calculated between subsequent images [START_REF] Michael | Introduction[END_REF]) exceeded 0.02 pixels were excluded before the point-by-point comparison with the numerical model [START_REF] Kok | Femoral strength and strains in sideways fall: Validation of finite element models against bilateral strain measurements[END_REF].

The following data comparison was developed according to published procedures on femur [START_REF] Kok | Femoral strength and strains in sideways fall: Validation of finite element models against bilateral strain measurements[END_REF][START_REF] Katz | New insights on the proximal femur biomechanics using Digital Image Correlation[END_REF]. A sphere, with a diameter of 1.0 mm, was calculated relative to the barycenter of each element of the analyzed surface meshed. All DIC strain measurement points within the sphere were averaged, and the obtained value was used for comparison to the FE element strain. In addition to a visual analysis, the minor and major principal strains were compared by performing a linear regression analysis and calculating the coefficient of determination (R²), the root mean squared error of estimate (RMSE), and the maximum error (Max err; maximal absolute difference between experimental and numerical values). Dispersion was evaluated using Bland-Altman plots [46]. Analysis was performed for each density-elasticity relationship considered in the current study.

To assess the influence of the ±2° uncertainty during model alignment, the average value of the Von Mises strain was also computed and compared to the experimental values.

Failure Load Comparison

The experimental reaction load curves over time were recorded to retrieve the maximum load of each fractured specimen during the impact. The obtained failure loads were compared with the simulations for each density-elasticity relationship and failure criteria. Due to the limited number of samples, non-parametric statistical tests were performed for this analysis by calculating Spearman's rank correlations ρ, RMSE, and Max err. The dispersion between these values was evaluated using Bland-Altman plots [46].

Fracture Prediction

Bone fracture was assessed by using the high-speed recordings, but also by a radiologist who interpreted radiographs after the impact. The medians, min, and max of numerical failure loads were calculated. The difference between the fractured and non-fractured bone maximum loads was evaluated using an unpaired Mann-Whitney U test. The small number of samples in each gender did not allow for a separation between males and females for the fracture discrimination analysis, even if significant differences were found between sexes.

A sensitivity/specificity analysis was performed by plotting Receiver Operating Characteristics curves (ROC curves) with their area under the curve (AUC) to determine the optimal threshold for each outcome. This optimal threshold was defined as the furthest point of the curve from the diagonal corresponding to the maximum of the sensitivity + specificity index. For each optimal detection threshold, test accuracy was calculated as the ratio of true positives (TP) and true negatives (TN), divided by the number of samples (

𝑇𝑃+𝑇𝑁 𝑛

). An analysis of variance was carried out on the AUCs compared to 0.5 in order to test whether each diagnosis test was more efficient than a simple random rule [START_REF] Hanley | A method of comparing the areas under receiver operating characteristic curves derived from the same cases[END_REF] and to evaluate the statistical differences between AUCs from the different density-elasticity relationships and failure criteria [START_REF] Delong | Comparing the Areas under Two or More Correlated Receiver Operating Characteristic Curves: A Nonparametric Approach[END_REF]. All tests were carried out using R software (R Foundation for Statistical Computing, Austria [START_REF] Rr Core | R: A language and environment for statistical computing[END_REF]).

Results

Experiments

Among the 30 radii, the experiments failed for two specimens due to misalignments of the potted bones, and the fracture was uncertain for an additional specimen. Finally, 27 bones were tested successfully. As described in the previous study [START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF]

Strain Comparison

The minor and major principal strains corresponding to the maximum experimental load were evaluated for each specimen and compared with the predicted strains of the finite element simulation. The sixteen non-fractured bones were considered for the first model evaluation. FE predicted major and minor principal strains had a strong correlation with the experimental measured strains (Figure 4) with R² = 0.76 for Eq. ( 1), R² = 0.78 for Eq. ( 2) and R² = 0.77 for Eq. ( 3). The slope of the linear regression was 1.06 for Eqs. ( 2)-(3) and 2.77 for Eq. ( 1). Moreover, Bland-Altman plots showed a bias of 0.27%ε for Eq. ( 1) with a ±2×SD of 1.73%ε and a bias of almost zero for Eqs. ( 2) and (3) with a ±2×SD of ∼0.40%ε.

In term of Von Mises strain average, the finite element predicted that strains vary as a function of Eqs. ( 1)-(3) (Figure 5). No significant difference was found between Eq. ( 2) and

(3) (p > 0.05) predicting more accurately the experimental values for non-fractured bones than Eq. ( 1), which overestimates the results. The ±2° uncertainty showed a variation of ±0.08% for Eqs. ( 2) and (3) and ±0.18% for Eq. ( 1) in the von Mises strain values. Correlation coefficients for experimentally measured strains versus predicted strains ranging from ρ = [0.68-0.71] (p < 0.01) for Eq. ( 1), to ρ = [0.75-0.82] (p < 0.001) for Eq. ( 2) and ρ = [0.78-0.79] (p < 0.001) for Eq. ( 3) (Table 2). The RMSE and Max Error were smallest for Eqs. ( 2) and ( 3)

and largest for Eq. ( 1) (Table 2). Qualitatively, similar results were found in the von Mises spatial strain distribution for the three equations with the min and max values of each field located at the same area. Furthermore, patterns present common characteristics with the experimental fields (Figure 6).

Failure Load Comparison

The eleven fractured bones were considered for the second model evaluation. The ±2° uncertainty was no longer considered with regard to the few effects on the surface von Mises strain for Eqs. ( 2) and (3). The numerical failure load varied as a function of Eqs. ( 1)-(3), failure criterion, and failure volumes (Figure 7). For the von Mises stress and Mohr-Coulomb criterion, where ultimate strength relied on hydroxyapatite equivalent density, computed failure loads increased with contiguous failed volumes. The lowest RMSEs were found for a volume of 450 mm 3 failed elements for these two criteria (ranging from 895 to 1181 N, corresponding to an error of 37 to 49% of the mean failure force measured experimentally).

Max error was smallest for Eq. ( 1) and largest for Eq. (2) (Supplementary Material).

Regarding the effective strain criterion based on Pistoia et al. [START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF], a fixed threshold of 0.7% strain leads to lower error for Eqs. 2 and 3. With this given criterion, the lowest RSME (558 N), max relative error (999 N), and average differences (-51 N) were found for a contiguous failed volume of 250 mm 3 and Eq. 2 corresponding to a mean error (RMSE%) of 23%. The Effective strain criterion, with a contiguous failed volume of 250 mm 3 and the densityelasticity equation 2, was considered in the rest of the study.

The correlation coefficient for experimentally measured failure load versus predicted failure load was ρ = 0.74 (𝑝 = 0.01). Data was compared with the literature to see if the results were part of a larger population. The radius failure loads of three different validated FE models with continuous meshing were extracted with their respective experimental data [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF][START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF][START_REF] Varga | Validation of an HR-pQCTbased homogenized finite element approach using mechanical testing of ultra-distal radius sections[END_REF] and plotted in Figure 8. These studies have shown good agreement between quasistatic experimental tests and FE analysis (R² = 0.87 [START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF], R² = 0.93 [START_REF] Varga | Validation of an HR-pQCTbased homogenized finite element approach using mechanical testing of ultra-distal radius sections[END_REF] and R² = 0.95 [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF]).

For an easier comparison with the literature results, the coefficient of determination was calculated for the current study (R² = 0.73, p < 0.05). Computed failure loads of the current study were in good agreement with the literature and can confirm good results for a larger population (Figure 8).

Fracture Prediction

A significant difference was found between experimental fractured and non-fractured bones maximum loads (p < 0.05). Median, min, and max are summarized in Table 3. A sensitivity/specificity study was performed by plotting ROC curves (Figure 9). The area under the curve (0.76) was significantly different from a random rule. An optimal threshold of detection was found with a sensitivity of 82% and specificity of 81%. This threshold allowed an accuracy of 82% for radius fracture prediction by subject-specific FE analysis in the case of a forward fall.

Discussion

This study investigated whether radius hFE modeling based on clinical imaging on a large distal area could discriminate ex vivo fractured bones from non-fractured bones. The model development required a sensitivity analysis on parameters such as the relation between the image grey level and elastic properties as well as failure criterion. The credibility of the model was assessed by comparing experimental and numerical surface strains and failure loads.

This is the first study to use enlarged distal radius hFE modelling to discriminate fractured from non-fractured bones.

Strain values in specific areas were similar between numerical simulations and experimental analysis on non-fractured bones for Eqs. ( 2) and ( 3) (R² = 0.77, p < 0.001, quadratic error of 21%ε with a slope of 1.06 for the principal strains. Spearman's ρ between 0.75 and 0.82, p < 0.0001, quadratic error between 0.14 and 0.19%ε with a slope of ∼1.0 for the average Von Mises strain). These results illustrate that experimental surface strains (induced by an impact at 2 m/s) in specific area in elastic response can be correctly predicted by quasi-static simulations. A point-by-point comparison (Figure 4) is given for the first time on the radius.

Such a comparison was proposed on the femur [START_REF] Kok | Femoral strength and strains in sideways fall: Validation of finite element models against bilateral strain measurements[END_REF][START_REF] Katz | New insights on the proximal femur biomechanics using Digital Image Correlation[END_REF] and even if a direct comparison cannot be done, the current results are in the same range (e.g., slope of the regression and mean values of the difference between simulated and experimental strains).

Literature showed that bone mechanical behavior is strain rate dependent. However, it is not a standard practice to use strain as a trustworthy parameter and the authors rather investigated the strain rate influence on stress or failure loads [START_REF] Jazinizadeh | Comparing the fracture limits of the proximal femur under impact and quasi-static conditions in simulation of a sideways fall[END_REF]. The authors have investigated the mechanical properties of cancellous bone (including strain) from quasi-static to dynamic strain rates with an emphasis on the intermediate regime [START_REF] Prot | The behavior of cancellous bone from quasistatic to dynamic strain rates with emphasis on the intermediate regime[END_REF]. The authors showed that within different regimes with a distinct strain rate from 0.001/s to 600/s, the yield stress and apparent modulus were strain rate dependent (significantly higher with higher strain rate). However, no significant difference was found for yield strain, suggesting that strain in elastic region is not greatly influenced by strain rate. The same results were found for the cortical bone [START_REF] Katzenberger | Effects of sex, age, and two loading rates on the tensile material properties of human rib cortical bone[END_REF]. These results suggest that strain should be similar between quasi-static and dynamic loadings.

The strong variation found for Eq. ( 1) illustrated the sensitivity of the model to the densityelasticity relationship. Because many of these equations exist in the literature [START_REF] Knowles | Quantitative Computed Tomography (QCT) derived Bone Mineral Density (BMD) in finite element studies: a review of the literature[END_REF], particular care should be made when developing a subject-specific finite element model. Some studies have investigated the influence of the density-elasticity relationship on the accuracy of specimen-specific finite element predicted strains [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Schileo | Subject-specific finite element models can accurately predict strain levels in long bones[END_REF][START_REF] Austman | The effect of the densitymodulus relationship selected to apply material properties in a finite element model of long bone[END_REF]. Edwards and Troy [START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF] compared three density-elasticity relationships for the human radius under quasi-static axial loading and found that Eq. ( 1) of the present study [START_REF] Keller | Predicting the compressive mechanical behavior of bone[END_REF] showed the best agreement between numerical and experimental surface strains. Nonetheless, this relationship provided the poorest agreement between experimental and predicted surface strains in the present study.

The most notable differences with the present study occurred between experiments. When Edwards and Troy loaded the radius axially through carpal bone and ligament at a low strain rate, the present study impacted the bones through a simplification of the radio-carpal joint at a high strain rate. This confirms that bone response is heavily dependent on boundary conditions such as the strain rate [START_REF] Carter | The compressive behavior of bone as a two-phase porous structure[END_REF][START_REF] Hansen | The Effect of Strain Rate on the Mechanical Properties of Human Cortical Bone[END_REF] or loading transfer [START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF]. Schileo et al. [START_REF] Schileo | Subject-specific finite element models can accurately predict strain levels in long bones[END_REF] also compared three density-elasticity relationships, although in human femurs, under several loading orientations and observed that Eq. ( 2) of the present study produced the best results.

rate. Similar studies show that with a higher displacement rate the fracture load should increase [START_REF] Jazinizadeh | Comparing the fracture limits of the proximal femur under impact and quasi-static conditions in simulation of a sideways fall[END_REF][START_REF] Courtney | Effects of loading rate on strength of the proximal femur[END_REF]. An increased fracture load could also be expected for the radius. Thus, the fracture loads of the current study should have been in the upper range of the Figure 8, including studies performed in quasi-static conditions. However, other parameters will affect the fracture loads, such as sex, age, size of the bone, bone mineral density, and orientation of the loading. In Figure 8, the results from Hosseini et al. 2017 [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF] present the highest values.

The failure loads were obtained on radius sections (pure compression) from younger subjects (77.5±9 y.o.) compared to other studies (81±10 y.o. for [START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF], 81±14 y.o. for [START_REF] Varga | Validation of an HR-pQCTbased homogenized finite element approach using mechanical testing of ultra-distal radius sections[END_REF], and 79.7±13 y.o. for the current study) and with a larger proportion of males. This may explain the larger values for that study. Moreover, the results from these studies were obtained on distal radius loaded in the direction of the diaphysis axis ("pure compression"). On the contrary, the current study considered an angle between the loading direction and the diaphysis axis. This loading configuration reduces the failure load compared to the loading in "pure compression". This important difference is probably the major reason why the failure loads in the current study are in the lower range of the failure loads reported in Figure 8, despite the impact loading. The majority of the fractured bones of the current study was considered osteoporotic by ultra-distal radius aBMD (DXA) [START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF]. This information was not available in the other studies. These low bone mineral densities might also explain the position of the current specimens compared to the literature. This comparison therefore highlighted that a sample size of 11 specimens for model validation is limited and this kind of analysis needs a larger population to confirm accuracy. However, the results suggested that the quasi-static hFE models of the distal radius can predict failure load obtained after an impact on the joint surface.

Computed failure loads on the current study computed by hFE analysis were associated with incident fracture. Predictions of HR-pQCT µFE and DXA measurement on the same pool of experimental data were previously established [START_REF] Revel | Influence of Loading Conditions in Finite Element Analysis Assessed by HR-pQCT on Ex Vivo Fracture Prediction[END_REF]. Derived µFE parameters showed an AUC of 0.795 with 82% accuracy for fracture detection of the best threshold (Sensitivity: 73%, Specificity: 88%). The DXA aBMD of the distal radius presented statistical similar results with an AUC of 0.773 and 85% accuracy (Sensitivity: 82%, Specificity: 88%) [ by including a second power term in the density-elasticity relationship [START_REF] Carter | The compressive behavior of bone as a two-phase porous structure[END_REF], could be also done to best fit the experimental conditions. Our group showed in a comparative study that the cortical bones from the radius and the femur present the same toughness under a strain rate representative of a fall [START_REF] Gauthier | Strain rate influence on human cortical bone toughness: A comparative study of four paired anatomical sites[END_REF]. This suggests that similar behavior for both sites is expected in the case of a fall-loading condition.

There are some limitations in the current study. First, the small sample consisted of 15 men and 15 women. Most of the non-fractured bones came from male donors, while the majority of fractured bones came from female donors [START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF], which may induce a bias during the sensitivity analysis. Moreover, the small number of samples did not allow for the possibility of full proper model validation, which requires a larger number of specimens. This problem was partially solved using previous published studies found in the literature, but different modelling techniques or experimental conditions remain present in these experiments.

Analyses of a larger number of single-sex samples would eliminate this bias. Other limitations are related to the use of cadaveric samples. Like most ex vivo studies on bones, the specimens were frozen, which is considered the best conservation mode. Moreover, the articular joint is simplified, and the load was evenly distributed on the articular surface of the radius while carpal bones have an effect on radius bone strength [START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF]. These simplifications differ from the in vivo phenomenon occurring during a fall in which soft tissue, ligament, and carpal bones are present and transfer the load to the radius.

Conclusions

Quasi-static hFE simulations could predict the elastic deformations of an external surface of the radius during dynamic experiments of forward falls on the forearm and depend on density-elasticity relationships. Failure criterion adapted from Pistoia et al. [START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF] was added to the linear elastic constitutive law to predict distal radius failure load. The originality of the current study is related to the evaluation of a homogenized finite element model to discriminate the two groups of bones (fractured and non-fractured) obtained through a unique loading condition. These results suggest that a credible FE modeling method is a suitable technique for predicting radius fractures occurring after a forward fall. Figure 6: Qualitative comparison of strain distributions for non-fractured specimens between A) experimental and computational model using B: Eq. ( 1), C; Eq. ( 2), and D: Eq. ( 3). Notice that the scale of colormap for Eq [START_REF] Khosla | A Crisis in the Treatment of Osteoporosis[END_REF] in B is twice as wide than the other to allow visualization of the strain field. Current study and data from literature [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF][START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF][START_REF] Varga | Validation of an HR-pQCTbased homogenized finite element approach using mechanical testing of ultra-distal radius sections[END_REF]. SUPPLEMENTARY Table 4 : For the 16 non fractured specimens: Baseline (mm), distance between cameras. Stereo angle (°) between cameras. Distance to specimen (mm) between cameras and the surface of the radii. This distance was not specifically measured. AOI Hight (mm and px), the long dimension of the analyzed region. Avg projection error (px) after calibration. Facet size or Subset (px) was chosen based on the average confidence interval of the correlation tracking function (sigma estimate ≈0.01 pixels).

Step (px, mm), refer to the overlap, was roughly a quarter of the subset in agreement with manufacturer guidelines. Filter size (px, mm) was chosen constants for all specimens. The noise trough time (%def) describes the spread of strain between consecutive images before the impact. Sigma Mean (px), the average error (sigma) of the analyzed images. Step (overlap, px)

Step

Figure 1 :

 1 Figure 1: Experimental setup with: (A) 3D scan of the experimental articular surface, (B) Resulting 3D surface of the bone, (C) Experimental setup, and (D) Capture from the front camera with its surface strain analysis.

Figure 2 :

 2 Figure 2: For the three density-elasticity relationships investigated, representation: (Left) Plot of Young's modulus against 𝝆 𝑯𝑨 . (Right) Longitudinal cross section illustrating material properties distribution.

Figure 3 :

 3 Figure 3: Top -Experiment setup from Zapata et al. (2017)[START_REF] Zapata | An ex vivo experiment to reproduce a forward fall leading to fractured and non-fractured radii[END_REF]. Bottom -Finite element model orientation and loading conditions.

Figure 4 :

 4 Figure 4: Left -Predicted vs. measured principal strains for each non-fractured specimen. Right -Bland Altman plots of predicted versus measured principal strains for the three density-elasticity relationships investigated.

Figure 5 :

 5 Figure 5: Left -Predicted (considering an uncertainty of ±2° for the radius orientation in regards to the load direction) vs. measured mean von Mises strains for each of non-fractured specimen. Right -Bland Altman plots of predicted versus measured von Mises strains.

Figure 7 :Figure 8 :

 78 Figure 7: Prediction error (RMSE) of the numerical model for the three different failure criteria, three elasticdensity relationships (Eqs. (1)-(3)) and for the four contiguous failed volumes considered.(Bottom Right) The Bland-Altman plot of the highest predicted failure loads obtained with effective strain criterion, 250 mm 3 failed volume and Eq. (2). The average difference is -51N Figure 8: Scatter plots between experimental and computed failure loads. (Left) Only the current study; (Right)

Figure 9 :

 9 Figure 9: ROC curves (Sensitivity against 1 -Specificity) of the computed failure load to predict experimental bone fracture. The dotted line represents the best threshold detection. *Significantly different from a random rule, alpha level = 0.05

Figure 1 Figure 2 Figure 3 Figure 5 Figure 6

 12356 Figure 1

Figure 7 Figure 8

 78 Figure 7

Figure 9 .

 9 Figure 9.

  ). Effective strain based on Pistoia et al. was calculated from the strain-energy density (U) and the

	Young's modulus E of each element: 𝜀 𝑒𝑓𝑓 = √ 2𝑈 𝐸	with a compressive yield strain (𝜀 𝑦𝑐 ) set to
	0.7% [40]. As inexistant on radius, the maximal von Mises stress and Mohr-Coulomb theories
	were based on compressive (𝜎 𝑦𝑐 ) and tensile (𝜎 𝑦𝑡 ) yield strengths calculated from
	hydroxyapatite equivalent density using pre-established equations on femoral cortical bone:
	𝜎 yc = 𝜎 yt = 137 × (1.034 × 𝜌 𝐻𝐴 + 0.0546) 1.88	for	0 < 𝜌 𝐻𝐴 ≤ 0.35	[38,41]
	𝜎 y𝑐 = 96.9 × 𝜌 𝐻𝐴 -5.95 and	𝜎 y𝑡 = 64 × 𝜌 𝐻𝐴 -21.82 for	𝜌 𝐻𝐴 > 0.35 [36,41]

  [START_REF] Revel | Influence of Loading Conditions in Finite Element Analysis Assessed by HR-pQCT on Ex Vivo Fracture Prediction[END_REF]]. The current model showed similar fracture prediction and was comparable to the technique used in the clinical research. If DXA has little room for progress, as well as µFE, limited by acquisition in vivo and restricted area, the hFE analysis with a clinical scan at a lower resolution has a large field of improvement. For further development, an enlarged distal radius FE model would allow more complex boundary conditions compared to a restricted area (e.g., different loading cases occurring in vivo). With this type of modelling, a full wrist FE model is conceivable with the presence of carpal bone, cartilage, and tendons which can apply stress to the bone. Moreover, we determined that the numerical response depends on densityelasticity relationships; no such type of relationship was developed specifically for the distal radius. Its development in fall conditions should help improve FE model accuracy. Inclusion of strain-rate-dependent behavior in our model, which could be done for Young's modulus

Acknowledgements

The authors would like to thank Karine Bruyère-Garnier for her help for the use of Vic3D.

On the human ulna, Austman et al. [START_REF] Austman | The effect of the densitymodulus relationship selected to apply material properties in a finite element model of long bone[END_REF] tested six different density-elasticity relationships under a simplified cantilever bending scenario and found the best agreement between numerical and experimental data for the two relationships from Carter and Hayes [START_REF] Carter | The compressive behavior of bone as a two-phase porous structure[END_REF], which are not used in this study, and the pooled bone site relationship described by Morgan et al. [START_REF] Morgan | Trabecular bone modulus-density relationships depend on anatomic site[END_REF] with Eq. ( 2). The differences in accuracy between studies confirmed the dependance of anatomical site [START_REF] Morgan | Trabecular bone modulus-density relationships depend on anatomic site[END_REF] and boundary conditions [START_REF] Edwards | Simulating distal radius fracture strength using biomechanical tests: a modeling study examining the influence of boundary conditions[END_REF][START_REF] Troy | Off-axis loads cause failure of the distal radius at lower magnitudes than axial loads: a finite element analysis[END_REF][START_REF] Carter | The compressive behavior of bone as a two-phase porous structure[END_REF] of the density-elasticity relationships on model accuracy.

The current study suggested that the effective strain criterion based on Pistoia et al. [START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF] and currently used in µFE HR-pQCT analysis produced the best agreement (ρ = 0.74) with the lowest amount of error (23%). But the small number of samples did not allow observation of whether the FE method of the current study could be applied to an enlarged population.

Similar in vitro studies showed agreement with the squared Pearson's R² ranging from 0.72 to 0.95 [START_REF] Hosseini | Fast estimation of Colles' fracture load of the distal section of the radius by homogenized finite element analysis based on HR-pQCT[END_REF][START_REF] Arias-Moreno | Validation of distal radius failure load predictions by homogenized-and micro-finite element analyses based on second-generation high-resolution peripheral quantitative CT images[END_REF][START_REF] Van Rietbergen | A survey of micro-finite element analysis for clinical assessment of bone strength: the first decade[END_REF][START_REF] Varga | Validation of an anatomy specific finite element model of Colles' fracture[END_REF][START_REF] Edwards | Finite element prediction of surface strain and fracture strength at the distal radius[END_REF][START_REF] Varga | Validation of an HR-pQCTbased homogenized finite element approach using mechanical testing of ultra-distal radius sections[END_REF][START_REF] Burkhart | Development and validation of a distal radius finite element model to simulate impact loading indicative of a forward fall[END_REF], depending on the number of specimens, constitutive law, loading scenario(s), model meshing technique (voxel vs. geometry based), and failure criterion.

When available and with similar modelling approaches, failure loads from the literature were extracted and compared with the results of the current study. These data compared quasistatic numerical failure load with low strain rate experimental loading until fracture. Despite strain rate and impact influencing the mechanical behavior of the bone [START_REF] Pithioux | Comparison of compact bone failure under two different loading rates: experimental and modelling approaches[END_REF][START_REF] Enns-Bray | Strain rate dependency of bovine trabecular bone under impact loading at sideways fall velocity[END_REF][START_REF] Jazinizadeh | Comparing the fracture limits of the proximal femur under impact and quasi-static conditions in simulation of a sideways fall[END_REF][START_REF] Gauthier | Strain rate influence on human cortical bone toughness: A comparative study of four paired anatomical sites[END_REF]], the fracture loads of the current study were added to the scatter plots of the literature to see their position in an enlarged population and show the credibility of the model before looking at fractures discrimination. Computed failure loads of the current study were in good agreement with the literature even if values were in the lower side of the range. However, conflicting results have been reported on femurs regarding the effect of the displacement

Conflict of interest

There is no conflict of interest for any of the authors.

Table List

Table 1: The failure criterion investigated with their respective equations. 𝜺 𝒆𝒇𝒇 is the effective strain estimated from the strain energy density [START_REF] Pistoia | Estimation of distal radius failure load with micro-finite element analysis models based on three-dimensional peripheral quantitative computed tomography images[END_REF]; 𝛆 𝒚𝒄 the compressive failure strain; 𝝈 𝒗𝒎 the von Mises stress; 𝝈 𝟏 , 𝝈 𝟐 , and 𝝈 𝟑 are the principal stresses (𝝈 𝟏 > 𝝈 𝟐 > 𝝈 𝟑 ); and 𝝈 𝒚𝒄 and 𝛔 𝒚𝒕 the compressive and tensile failure stresses.

Criterions Equations

Effective strain

Maximal von Mises stress

Mohr-Coulomb SUPPLEMENTARY Figure 11 : Mesh convergence study on a sub-sample of 6 specimens ranging from 50,000 to 600,000 elements by increments of 50,000. The mean equivalent von Mises strain in a region of interest of the anterior face of the radius was evaluated.

SUPPLEMENARY