Cluster Analysis of Finite Element Analysis and Bone Microarchitectural Parameters Identifies Phenotypes with High Fracture Risk

High-resolution peripheral quantitative computed tomography (HRpQCT) is increasingly used for exploring associations between bone microarchitectural and finite element analysis (FEA) parameters and fracture. We hypothesised that combining bone microarchitectural parameters, geometry, BMD and FEA estimates of bone strength from HRpQCT may improve discrimination of fragility fractures. The analysis sample comprised of 359 participants (aged 72–81 years) from the Hertfordshire Cohort Study. Fracture history was determined by self-report and vertebral fracture assessment. Participants underwent HRpQCT scans of the distal radius and DXA scans of the proximal femur and lateral spine. Poisson regression with robust variance estimation was used to derive relative risks for the relationship between individual bone microarchitectural and FEA parameters and previous fracture. Cluster analysis of these parameters was then performed to identify phenotypes associated with fracture prevalence. Receiver operating characteristic analysis suggested that bone microarchitectural parameters improved fracture discrimination compared to aBMD alone, whereas further inclusion of FEA parameters resulted in minimal improvements. Cluster analysis (k-means) identified four clusters. The first had lower Young modulus, cortical thickness, cortical volumetric density and Von Mises stresses compared to the wider sample; fracture rates were only significantly greater among women (relative risk [95%CI] compared to lowest risk cluster: 2.55 [1.28, 5.07], p = 0.008). The second cluster in women had greater trabecular separation, lower trabecular volumetric density and lower trabecular load with an increase in fracture rate compared to lowest risk cluster (1.93 [0.98, 3.78], p = 0.057). These findings may help inform intervention strategies for the prevention and management of osteoporosis. Electronic supplementary material The online version of this article (10.1007/s00223-019-00564-7) contains supplementary material, which is available to authorized users.


Introduction
Although measurements of areal bone mineral density (aBMD) of the hip or spine is considered the gold-standard for assessing fracture risk [1], it is widely recognised that there are limitations to only using DXA to determine bone fragility. Stratifying for fracture risk using DXA-derived aBMD may fail to identify some individuals at higher risk of fracture; indeed, as many as half of hip fractures occur in individuals with aBMD values considered low risk for osteoporotic fracture [2]. Prediction tools including clinical risk factors and DXA-derived BMD, such as the Fracture Risk Assessment Tool (FRAX ® ), improve prediction of fractures [3] but still fail to identify many who go on to sustain a fracture [4]. Understanding better the bone phenotypes associated with bone fragility may improve fracture prediction. To address this, techniques such as high-resolution peripheral quantitative computed tomography (HRpQCT) have been developed to explore volumetric bone mineral density (vBMD) and bone microarchitecture in both trabecular and cortical bones. It is now also possible to perform finite element analyses on HRpQCT scans and previous research has related tibial and radial FEA parameters to increased risk of fragility fractures at all sites, independent of aBMD [5][6][7]. Several other studies have reported HRpQCT parameters to be associated with prior or future fracture, independent of aBMD [7][8][9][10][11].
Our previous work in the Hertfordshire Cohort Study (HCS) used cluster analyses to identify two similar phenotypes associated with increased risk of prevalent fracture [12]. One demonstrated a 'cortical deficiency' phenotype, characterised by lower cortical thickness and cortical volumetric density and in men only, higher total and trabecular area. Men in this cluster had a normal aBMD and, therefore, may not have been identified as high risk by conventional DXA-based risk stratification alone. Men and women in the other cluster showed a 'trabecular deficiency' phenotype, with lower trabecular volumetric density and number. Similar radial clusters were identified in the Global Longitudinal Study of Osteoporosis in Women (GLOW) [13].
The aim of the current study was therefore to determine the extent to which bone microarchitectural and FEA parameters improve fracture discrimination compared to using aBMD alone. We hypothesised that combining bone microarchitectural parameters, geometry, BMD and FEA estimates of bone strength from HRpQCT as a composite of bone strength may improve discrimination of fragility fractures. The secondary aim was to repeat the cluster analyses, to see whether FEA parameters changed the identified phenotypes previously associated with fracture.

The Hertfordshire Cohort Study
The HCS comprises 2997 men and women born in Hertfordshire from 1931-1939 and who still lived there in 1998-2004 when they completed a baseline home interview and research clinic for a detailed characterisation of their health; the study has been described in detail previously [14]. At the baseline home interview (1998)(1999)(2000)(2001)(2002)(2003)(2004), menopausal status (women only) and customary physical activity level (Dallosso questionnaire [15]) were ascertained by a nurse-administered questionnaire. Dietary calcium intake was determined using a food-frequency questionnaire [16]. Social class was coded from the 1990 OPCS Standard Occupational Classification (SOC90) unit group for occupation [17].
In 2011-2012, 570 East Hertfordshire participants were invited to take part in a further follow-up study; 376 agreed to participate [12]. Smoking status, alcohol consumption and whether participants had broken any bones since aged 45 years were ascertained by a nurse-administered questionnaire. Information was ascertained on whether participants had used bisphosphonates or hormone replacement therapy since HCS baseline.
The baseline HCS had ethical approval from the Hertfordshire and Bedfordshire Local Research Ethics Committee and the follow-up had ethical approval from the East and North Hertfordshire Ethical Committees. Investigations were conducted in accordance with the principles expressed in the Declaration of Helsinki.

Participant Assessments at the 2011-2012 Follow-up
On the day of scanning, height was measured (wall-mounted SECA stadiometer) along with weight (calibrated SECA 770 digital floor scales, SECA Ltd, Hamburg) and used to derive BMI (kg/m 2 ). Bilateral scans of the proximal femur were used for assessment of femoral neck aBMD (Lunar Prodigy Advance DXA scanner (GE Medical Systems)); the lowest value was used for analyses and diagnosis of osteoporosis or osteopenia. Morphometric vertebral fractures were diagnosed from a lateral spine view imaged using the same machine and graded based on the Genant semi-quantitative method of vertebral fracture assessment [18]. Participants with a vertebral fracture or a self-reported fracture since age 45 years were regarded as having had a previous fracture.
HRpQCT scans (XtremeCT, Scanco Medical AG, Switzerland) of the non-dominant distal radius were performed; dominant limbs were scanned if the non-dominant limb had previously fractured. In total, 110 parallel CT slices were obtained, representing a volume of bone 9 mm in axial length with a nominal resolution (voxel size) of 82 μm. The scan protocol was in accordance with manufacturer's guidelines and as described by Boutroy et al. [19]. Using the method of Pauchard et al. [20], eight scans were excluded due to excessive motion artefact (grade 5 scans); scans of quality 4 and above were included in the analysis. Manufacturer standard evaluation and cortical porosity scripts were used for image analysis [5,[21][22][23][24][25]. Cortical and trabecular densities described in this study were ascertained using HRpQCT and are volumetric (mg/cm 3 ).
A detailed description of the FEA analysis development for in vivo assessment of bone strength is in Boutroy et al. 2008 [25]. Briefly, the FE-solver (Image Processing Language) scripts were run to assess the biomechanical properties of the cortical and trabecular compartments and of the whole bone. The FEA models assume boundary conditions for an applied compressive load in the axial direction to the radius or tibia. From these the various stress, stiffness and failure load parameters were estimated.

Statistical Analysis
Data were described using summary statistics. Age, anthropometric and lifestyle characteristics were compared between individuals who did and did not have a previous fracture. Skewed bone microarchitectural and FEA parameters were transformed prior to standardising. Poisson regression with robust variance estimation was used to derive relative risks for the relationship between individual parameters and previous fracture. Unadjusted and fully adjusted relative risks, accounting for age, height, BMI, dietary calcium, physical activity, smoking history (ever vs never), alcohol consumption, social class, bisphosphonate use, time since menopause (women only) and hormone replacement therapy (women only), were estimated.
Risk of previous fracture was examined using sex-specific logistic regression models containing the following sets of predictors: femoral neck aBMD only; femoral neck aBMD and bone microarchitectural parameters that were significantly associated with fracture in fully adjusted analysis; FEA parameters that were significantly associated with fracture in sex-specific analysis as additional predictors. The following pairs of predictors were highly collinear and, therefore, were not included in the same model: trabecular and total area; bone stiffness and failure load; trabecular strain and von Mises stresses (trabecular). Performance of models was assessed using the areas under the receiver operating characteristic curves (AUC).
A cluster analysis of the complete set of bone microarchitectural and FEA parameters was performed using the k-means partitioning method. The number of clusters selected was based on the stability of the clustering, and on the potential for identifying contrasting phenotypes as in previous analyses [12,13]. The means and standard deviations (SD) of the standardised parameters, and fracture proportion were then determined for each cluster. Poisson regression with robust variance estimation was used to determine the likelihood of fracture in each cluster compared to the lowest risk cluster. Mean femoral neck aBMD in each cluster was compared to the cluster with the lowest fracture risk.
The analysis sample comprised of 359 participants with complete data for all radial bone microarchitectural and FEA parameters. Healthy participant effects were assessed by comparing HCS baseline participant characteristics between this analysis sample of 359 participants and the group of 2638 participants who attended the HCS baseline clinic but were not included in the analysis sample. All analyses were performed among men and women separately using Stata 15 (StataCorp. 2017. Stata Statistical Software: Release 15. College Station, TX: StataCorp LLC).

Participant Characteristics
The characteristics of the study population are presented in Table 1. Mean (SD) age of the 359 participants at the time of scan was 76.3 (2.6) years. Overall 45 (26.2%) men and 50 (31.6%) women had a previous fracture (vertebral or selfreported since age 45 years). Locations of fractures in this cohort have been described previously [12].
On average, women who had a previous fracture were older than those who did not (77.3 years vs. 76.1 years, p = 0.005); no significant differences were observed among men. Among men and women, the following characteristics did not differ significantly between those who did and did not have a fracture: height, weight, dietary calcium, physical activity, smoking status and alcohol consumption (data not shown).

Assessing Healthy Participant Effects in the Analysis Sample
Compared to the 2638 participants who attended the HCS baseline clinic but were not included in the analysis sample, both men and women in the analysis sample had higher self-reported physical activity at baseline (p < 0.006); only men in the analysis sample were more likely to have never smoked at baseline (p = 0.003) and only women were less likely to have high alcohol consumption (p = 0.047). The proportion who were of manual social class (classes IIIM, IV and V) did not differ significantly (p > 0.05) between the two groups.

Relationships Between Bone Microarchitectural and FEA Parameters and Previous Fracture
The associations between individual bone microarchitectural and FEA parameters and previous fracture among men and women are presented in Table 2. Among men, higher total and trabecular area and lower cortical thickness were each associated with increased risk of previous fracture in unadjusted and adjusted analyses (p < 0.03); no FEA parameters were associated with fracture risk (p > 0.05). Among women, lower values of the following parameters were associated with increased risk of fracture in unadjusted and adjusted analyses (p < 0.04): cortical area and porosity; trabecular density, thickness, Von Mises stresses and strain; bone stiffness and failure load and Young Modulus.

Comparison of Fracture Prediction Models Using Receiver Operating Characteristic Analysis
Areas under the receiver operating characteristic curves (AUC) for different models used to examine risk of fracture are presented in Table 3 Table 3 are stated in Online Appendix 1.

Cluster Analysis of Bone Microarchitectural and FEA Radius Parameters
Four clusters were obtained among men and women. Summary statistics of the standardised bone microarchitectural and FEA parameters, femoral neck aBMD and fracture prevalence according to the different clusters are illustrated in Tables 4 and 5. Figure 1 displays means of the standardised parameters and fracture prevalence for Clusters 1 and 2, the highest risk clusters. Unless otherwise indicated, statements about the bone microarchitectural and FEA parameters refer to sex-specific means which differed by more than one SD compared to the analysis sample; comparisons of aBMD and fracture prevalence are in relation to Cluster 4, the cluster with lowest fracture risk. Men in Cluster 4 had higher cortical thickness and Young modulus and lower trabecular area; women had greater Young modulus and stiffness (all differences in means > 0.9 SD). In Cluster 1, lower means of the following parameters were found among men and women: cortical thickness, density, Young modulus and Von Mises stresses. Among men only, cortical strain was also lower and among women only, Cluster 1 was associated with lower cortical area, trabecular density and thickness, bone stiffness and failure load. Fracture risk was only significantly greater among women (relative risk [95% CI]: 2.55 [1.28, 5.07], p = 0.008).
In Cluster 2, a trabecular deficient phenotype among women only was found, with higher trabecular separation, lower trabecular density and lower proportion of load applied to both distal and proximal trabecular bones (differences in means > 0.9 SD). This phenotype was associated with increased fracture risk in women (1.93 [0.98, 3.78], p = 0.057).
In Cluster 3, bone stiffness and failure load were higher among men only (differences in means ≥ 0.9 SD), though aBMD was not different to Cluster 4 (p = 0.108). No significant differences regarding fracture risk were observed.

Discussion
This study assessed the additional utility of adding HRpQCT, BMD, microarchitectural and FEA parameters to enhance fracture discrimination in a cohort of elderly men and women. Whilst individual measures of bone microarchitecture and FEA discriminated fracture cases versus nonfracture groups, and were selected in cluster analyses, little benefit of adding the FEA parameters was found over and above femoral neck aBMD and bone microarchitecture in terms of fracture discrimination.
In terms of individual bone microarchitectural parameters, higher total and trabecular area and lower cortical thickness were associated with increased fracture risk in unadjusted and adjusted analyses among men; corresponding parameters among women were lower cortical area and porosity, trabecular density and thickness. The cortical porosity result in women seems a little counterintuitive though is consistent with previous work in this cohort [12], and in GLOW [13]. The observation may be due to the way the cortical porosity analysis script defines a pore [26] which is based on how many neighbouring voxels have a similar low attenuation value. If an individual has less cortical bone and thinner cortices (as shown in our current analyses), it may be that fewer pores meet these criteria, which would be reflected in a lower % porosity. Secondly, the bone of participants with fracture may have a lower turnover, and thus repair rate, which may also result in fewer pores and an increase in fracture risk. In men, no FEA parameters were associated with fracture risk, but in women, lower bone stiffness, failure load, Young modulus and trabecular stress and strain were all associated with higher fracture risk in both unadjusted and adjusted analyses. It appears that men at higher risk of fracture had larger bones with a thinner cortex and deterioration of trabecular bone, which was not translated into reduced bone strength as measured by linear FEA. Linear FEA has been criticised as it only evaluates stresses and strains placed on the bone in one direction, representing a linear compressive force on the bone, whereas non-linear FEA may provide additional information on fracture risk [27]. Perhaps this phenotype in men may produce vulnerability to forces placed on bones from other directions, including bending, which would not be identified by linear FEA. If this were the case, it may be expected that these men may be at higher risk of hip or radial fractures, frequently sustained on falling, compared with vertebral, classically compressive, fractures. During ageing, men also compensate for bone loss with greater periosteal formation to a greater extent than women [28], which may protect against compressive forces better in males. Conversely, the impaired linear bone strength parameters seen in women may suggest these women are more prone to compressive fractures. Vertebral fractures were the most common fracture among both males and females, although due to the small number of fractures in our cohort, comparison of fracture types in different clusters was not possible. It would be interesting to investigate this further, using non-linear FEA, in a larger sample with larger numbers of fractures.
The ROC curve analyses showed the benefit of adding relevant bone microarchitectural parameters into DXA-derived femoral neck aBMD assessment of fracture risk. This benefit was not statistically significant, but this may be influenced by our small sample size. Despite the association of FEA parameters with fracture risk in females, ROC curve analysis revealed FEA had little additional benefit in predicting fracture over models using bone microarchitectural parameters and aBMD. FEA is computationally intensive, and therefore time consuming and expensive. Our analyses suggest the additional information may not justify its use in predicting fractures in clinical practice. These findings are in contrast with previous studies which have found FEA to be a valuable predictor of prevalent and incident fracture, more important than other HRpQCT measures [5][6][7]29]. One study found a machine learning model incorporating HRpQCT measures could predict fractures better than aBMD [30]. Another study reported the additional benefit of a combination of a trabecular and cortical parameter or FEA failure load to femoral neck aBMD or FRAX-BMD in predicting fracture in postmenopausal women. However, failure load could be replaced with aBMD of the ultra-distal radius with no significant reduction in predictive worth [24]. This is similar to our findings, with the addition of FEA providing little additional benefit over bone microarchitectural parameters.
Our cluster analysis revealed similar clusters to previous findings [12] with Cluster 1 showing a predominantly 'cortical deficiency' phenotype based on bone microarchitectural parameters. In males, trabecular area was greater, perhaps reflecting a greater proportion of trabecular bone due to a reduction in cortical bone. In females, trabecular density and thickness were also lower, suggestive of more generalised deterioration of bone structure. FEA in this cluster showed greater percentage of load on trabecular bone, lower Young modulus and lower cortical stresses in both males and females and lower bone stiffness and bone failure load in females. Compared to the lowest risk cluster, this cluster also had lower mean aBMD, and was the only cluster with significantly greater fracture risk in females.
The trabecular phenotype was less definitive than in our previous analysis but was indicated in Cluster 2 where females tended towards a 'trabecular deficiency' phenotype, and tended towards higher fracture risk, although this did not reach statistical significance. This is likely due to the additional FEA parameters included in the cluster analysis compared to the previous analysis which only included bone microarchitectural parameters. Bone microarchitectural and FEA parameters were similar to the wider sample, even though there was a significantly lower aBMD in both males and females compared to the lowest risk cluster.
These findings demonstrate the importance of HRpQCT parameters and identify cortical deterioration as a key contributor to fracture risk. As with our findings, previous studies have found both cortical and trabecular microarchitectural deterioration important for fracture risk, some finding cortical [31,32], and some finding trabecular [9, 19, 33] changes more important. One study found both lower trabecular and cortical volumetric density were independently associated with fracture incidence [24]. Our use of cluster analysis helps to elucidate different phenotypes within the population at risk of fractures, with different contributions of both trabecular and cortical parameters. The strengths of our study include basing our analyses on the HCS, a well characterised cohort where data were rigorously collected by an experienced multidisciplinary team. Furthermore, our analyses used both bone microarchitectural and FEA parameters to provide a comprehensive illustration of bone phenotypes.
This study has some limitations. Firstly, a healthy responder bias has been observed in HCS and examining participant characteristics according to inclusion status has revealed healthier lifestyles at baseline for participants included in the analysis sample compared to those who were not. However, our analyses were internal, so bias would only arise if the associations of interest differed systematically between those who were included in the analysis sample and those who were not; this seems unlikely. Secondly, temporal causation cannot be inferred as our study has a cross-sectional design. It may be that the differences in bone microstructure seen are secondary to remodelling in response to fracture, rather than properties of the bone which predispose to fracture, especially as we have only collected information about previous fractures. Thirdly, fracture status was missing for some participants, although this information was available for the vast majority (91.9%) of the analysis sample. Finally, the low numbers of reported fractures and a relatively small sample size, along with the lack of stability regarding cluster analysis algorithms in general, may limit the generalisability of findings. However, the similarity of the clusters observed to those in other analyses and their biological plausibility suggests that they are robust.
In conclusion, microarchitectural deterioration, bone geometry and, in women, FEA-derived bone strength contributed to an increased risk of previous fracture. Cluster analysis revealed a cortical and a trabecular deficiency phenotype, which both showed lower aBMD in men and women. Only women with the cortical deficiency phenotype had significantly increased risk of previous fractures. In this cohort, adding bone microarchitectural parameters to aBMD could better predict previous fracture, but further addition of FEA conferred little benefit.
Author Contribution CC (guarantor), EMD and KAW designed the study. LDW conducted the statistical analyses and drafted the methods and results. CS drafted the introduction and discussion. EMD and KAW are also investigators of the Hertfordshire Cohort Study. KAW and MHE provided guidance regarding the analysis strategy and implications of the findings. All authors made substantial contributions to the manuscript and approved the final version.

Compliance with Ethical Standards
Conflict of interest CC reports personal fees (outside the submitted work) from Alliance for Better Bone Health, Amgen, Eli Lilly, GSK, Medtronic, Merck, Novartis, Pfizer, Roche, Servier, Takeda and UCB. EMD reports personal fees (outside the submitted work) from Pfizer Fig. 1 Means of standardised radial parameters within Clusters 1 and 2. *Significantly higher fracture prevalence (p = 0.008) compared to Cluster 4 (lowest risk cluster). Ct cortical; Tb trabecular; Ct density cortical density; Tb density trabecular density; dist distal; prox proximal; mod modulus; Tb stress Von Mises stresses (trabecular); Ct stress Von Mises stresses (cortical) Healthcare and from the UCB Discussion panel. MHE received (outside the submitted work) an unrestricted project grant from Servier; conference attendance support from Eli Lilly, UCB and Pfizer and course attendance support from Chugai and Abbvie. LW, CS and KAW declare that they have no conflicts of interest.

Ethical Approval
The baseline Hertfordshire Cohort Study had ethical approval from the Hertfordshire and Bedfordshire Local Research Ethics Committee and the follow-up had ethical approval from the East and North Hertfordshire Ethical Committees.

Human and Animal Participants
All procedures performed in studies involving human participants were in accordance with the ethical standards of the institutional and/or national research committee and with the 1964 Helsinki Declaration and its later amendments or comparable ethical standards.
Informed Consent All participants gave signed consent to participate in the study and for their health records to be accessed in the future.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.