A method and device for determining a probability of fracture of a metastatic bone

Abstract

Metastases increase the risk of bone fracture. Consequently, clinicians need to know if the patient's metastatic bones can withstand the stress of daily activities, in particular femurs and vertebras. Accordingly, it is disclosed a method for determining a probability of fracture of a metastatic bone of a patient in daily activities, that is determined from a numerical simulation of load application performed on a Finite Element model of the bone, and from a determined uncertainty related to the building of the FE model and the numerical simulation. The uncertainty has been elaborated by studies of replicability and sensitivity of models. According to the disclosed method, various probabilities may be computed for different types of daily activities according to the average load represented by each activity on the bone. The method may be applied to metastatic bones or vertebras.

Description

[0001] A METHOD AND DEVICE FOR DETERMINING A PROBABILITY OF FRACTURE OF A METASTATIC BONE TECHNICAL FIELD The present disclosure relates to method and device for determining a risk of fracture of a metastatic bone of a patient, consecutive to a load being applied to the metastatic bone. This disclosure may in particular be applied to metastatic femurs or metastatic vertebrae. TECHNICAL BACKGROUND It is estimated that cancer affected about 19.3 million people in 2020. Around 80% of patients with malignant tumours develop bone metastases. Bone metastases are tumours arising from the migration of cancer cells from the primary cancer to a bone site. Bone metastases affect the mechanical strength of the bone and can, in extreme case, be the main cause of fracture during daily activities. Estimating the risk of fracture of a metastatic bone is essential for a medical practitioner in order to determine the best course of the cancer treatment. Indeed, if the risk of fracture during daily activities is high, the practitioner will opt for a prophylactic operation to avoid the complications of a fracture. However, this operation will require the patient to stop taking treatments that prevent cancer form progressing during the post-operative convalescence period. The risk of fracture must therefore be assessed as accurately as possible, so that surgery can be avoided until it is necessary. In this perspective, it is known from [Eggermont, 2022], a method applied to metastatic femurs, comprising building a Finite Element model from a CT-scan of the metastatic femur, performing numerical simulation using a Finite Element model until fracture of the femur. The strength of the femur is determined by the maximum total reaction force in Newtons, and the so-called Bone Strength (BOS) Score is computed by dividing the strength of the femur by the body weight of the patient. With reference to figure 1, the BOS score was compared to clinical data of patients having undergone, or not, a fracture, to determine an indicative threshold of a risk of fracture. The threshold was set to 7.5 times the body weight of the patient. As shown in figure 1, the BOS score enables to conclude that BOS scores greater than the threshold correspond to low fracture risk, but having a score lower than the threshold does not enable discriminating situations. In the publication by [Sas, 2022], the threshold for determining a risk of fracture was reevaluated to 8.9. The threshold thus appears different according to the considered population, and appears as insufficient for reliably evaluating a risk of fracture for a given patient. SUMMARY OF THE INVENTION The present disclosure aims at improving the situation. In particular, an aim of the present disclosure is to propose a method enabling an accurate assessment of a risk of fracture of a metastatic bone for a determined load value applied to the bone. Another aim of the present disclosure is to enable assessing a risk of fracture of a metastatic bone for a plurality of daily situations or activities of the patient, in order to provide a medical practitioner with relevant hindsight for choosing whether or not a surgery of the metastasis should be done. To this end, a computer-implementing method for determining a risk of fracture of a metastatic bone of a patient, consecutive to a load being applied to the metastatic bone, wherein the bone is either a femur or a vertebra is disclosed, the method comprising: - Obtaining a Computed Tomography (CT) scan of the bone, the CT scan comprising a plurality of voxels associated to respective density values, - Building, from the CT scan of the bone, a Finite Element Model, comprising: o segmenting the CT scan to extract a three-dimensional volume corresponding to the bone, o generating a three-dimensional meshed model of the bone, comprising a plurality of elements, from the three-dimensional volume, o attributing mechanical properties to the meshed-model, said mechanical properties being derived from the density values of the CT scan, - performing numerical simulation of load application on the FE model of the bone, according to pre-determined boundary conditions, to determine a simulated failure load value, and - determining a probability of fracture of the bone, for a determined load value applied to the bone, based on the simulated failure load value and a determined uncertainty related to the building of the FE model and the numerical simulation. In embodiments, the boundary conditions of the numerical simulation comprise parameters defining the position of load application, and a relative orientation of the bone with respect to the direction of the applied loads. In embodiments, the bone is a femur, and the boundary conditions of the numerical simulation comprise at least one of: - A length of a portion of femur shaft in the three-dimensional volume corresponding to the femur, - An orientation of the femur in a front plane, - A rotation angle of the femur about the femur shaft, - A position of application of the loads on the femoral head. In embodiments, the bone is a vertebra, and the boundary conditions of the numerical simulation comprise at least one of: - An orientation of the vertebra in a sagittal plane, - An orientation of the vertebra in a front plane, - A position of application of the loads on the vertebra. In embodiments, the uncertainty related to the building of the FE model results from derivation of mechanical properties from the density values of the CT scan and / or from the segmentation of the CT scan. In embodiments, the method comprises a preliminary step of determining an uncertainty value related to the building of the FE model and the numerical simulation, comprising: - For a plurality of ex-vivo bones: o Performing a CT-scan on the ex-vivo bone, building a Finite Element model of the ex-vivo bone, and implementing a numerical application of load application on the FE model of the ex-vivo bone until a predetermined failure criterion is fulfilled for a simulated failure load, the numerical simulation of load application being performed according to determined boundary conditions on the position of load application and the relative orientation of the bone with respect to the direction of the applied loads, o Performing load application on the ex-vivo bone according to the same boundary conditions as used in the numerical application, to determine an actual failure load, and o Computing the uncertainty value from the standard deviation between the simulated failure loads and the actual failure loads for the plurality of ex-vivo bones. In embodiments, determining a probability of fracture of the bone for a determined load value applied to the bone, from the simulated failure load value and the determined uncertainty value, comprises: - Determining, from the determined uncertainty value and from the simulated failure load value obtained from the numerical simulation performed on the bone, a density of probability of the failure load, and - Computing, for the determined load value applied to the bone, a probability of fracture as the integral of the density of probability up to the determined load value. In embodiments, the determined load value is a multiple of the body weight of the patient. In embodiments, the segmentation is performed by implementation of an automatic segmentation algorithm. In embodiments, the method comprises computing a probability of fracture of the bone, and associated uncertainty, for each of a plurality of determined load values corresponding to respective average load values applied to the bone in a plurality of situations or activities to which the patient may be subject, and displaying, for each situation or activity, the corresponding probability of fracture. In embodiments, the bone is a metastatic femur, and the method is applied to the metastatic femur of a patient and to the other femur of the patient, and further comprises displaying, for each femur of the patient, the probability of fracture of each femur for a same applied load value. According to another object, a computer-program product is disclosed, comprising code instructions for implementing the method according to the description above, when it is implemented by a computer. According to another object, a device for determining a risk of fracture of a metastatic bone of a patient, consecutive to a load being applied to the metastatic bone, wherein the bone is either a femur or a vertebra, the device comprising at least a computer and a memory, and being characterized in that it is configured for implementing the method according to the description above. In embodiments, the device further comprises a CT-scanner in communication with the computer, the CT-scanner being configured to acquire CT-scans and to communicate the acquired CT-scans to the computer. The method enables computing the probability of fracture of a metastatic bone, for a determined load value applied to the bone, based on the results of a numerical simulation performed on a Finite Element model of the bone, by taking into account the uncertainties related to the building of the Finite Element model and the boundary conditions of the numerical simulation. The method in particular enables computing the probability of fracture of a metastatic bone for a plurality of load values which may correspond to average load values exerted on the bone during a plurality of determined situations or daily activities such as: walking, sitting down, standing up, going down stairs, etc. BRIEF DESCRIPTION OF THE DRAWINGS Other features and advantages of the invention will be apparent from the following detailed description given by way of non-limiting example, with reference to the accompanying drawings, in which: - Figure 1, which has already been described, represents an example of Bone Strength Scores computed in the prior art document discussed above for a plurality of patients with metastatic femurs. - Figure 2a schematically represents the main steps of a method for determining a simulated failure load value of a bone, - Figure 2b schematically represents the main steps of a method according to an embodiment, - Figure 3 represents an exemplary device for determining a risk of fracture of a metastatic bone according to embodiments of the disclosure. - Figure 4 represents an example of non linear stress-strain curve that may be used for populating a finite element model with mechanical parameters ([Keyak, 2005]), - Figure 5a is an illustration of boundary conditions of load application on a femur, - Figure 5b is an illustration of boundary conditions of load application on a vertebra, - Figure 6a represents an example of density of probability of fracture of a femur and a probability of fracture at a load value equal to 3.5 times the body weight, corresponding to the load applied on the femur while descending stairs, - Figure 6b represents an example of cumulative fracture probability distribution function with load value. - Figure 7 represents an exemplary summary that can be output for a medical practitioner. DETAILED DESCRIPTION OF AT LEAST ONE EMBODIMENT A method for determining a risk of fracture of a metastatic bone of a patient will now be described. With reference to figure 3, this method is implemented by a device 1 comprising a computer 10, for instance processor, microprocessor, controller or microcontroller, and a memory 11 storing code instructions for implementing the method described below when they are executed by the computer. The memory may also store preliminarily determined parameters for implementing the method, and in particular the parameters for computing a probability of fracture of a bone from a simulated failure load value, as discussed in more details below. The memory may include a non-volatile electronic memory, an optical and / or magnetic drive, etc. In embodiments, the device 1 may also comprise a Computed Tomography scanner 20, configured for acquiring CT scans of bones, either in-vivo or ex-vivo. As commonly known, a CT scanner is a medical imaging device using X-rays and computational treatments to create cross-sectional images of a body. X-ray beams are emitted from different angles and pass through the body. Detectors located on the opposite side of the body with respect to the X-ray beam emitters measure the amount of radiation that have passed through. The amount of radiation absorbed by the body’s tissues varies depending on their density. The measured amounts of radiations are then processed to produce cross-sectional images of the body at different depths, which can be viewed individually or stacked together to create a 3D representation including a plurality of voxels, i.e. of three-dimensional pixels, of the area being examined. The voxels are associated to respective density values, which are inferred from the measured amounts of radiations having crossed the body. The Computed-Tomography scanner may be configured for transmitting acquired CT scans to the computer 10. The CT scanner and computer may be collocated, or the computer 10 may be located at a remote location relative to the CT scanner. In any case, the CT scanner may be configured to communicate with the computer using wired or wireless communications. With reference to figures 2a and 2b, are shown the main steps of a method for determining a risk of fracture of a metastatic bone of a patient according to embodiments of the disclosure. The metastatic bone may be a femur or a vertebra. The method may in particular be applied to an in-vivo metastatic bone, in order to provide a medical practitioner the assessment of a risk of fracture of the bone in daily activities of the patient. From said assessment, the practitioner may determine whether or not a surgery for removal of the metastases should be performed prior to continuing the cancer treatment. The method comprises a step 100 of the computer 10 obtaining a CT scan of the bone, the CT scan being a three-dimensional representation of the area of the body including the bone, comprising a plurality of voxels where each voxel is associated to a respective density value of the tissues contained within said area. The density value of the tissues that is assigned to each voxel may be computed from the intensity or grey level values of the voxels, according to a predetermined calibration. For instance, a calibration phantom may be positioned underneath the patient while acquiring the CT scan, where the density of the calibration phantom is known and enables establishing said calibration law between the grey level or intensity values of the voxels and density values. The CT scan may have been preliminarily acquired on a patient using a CT scanner 20, and stored in a memory for later processing, in which case step 100 comprises the computer 10 recovering said CT scan from the memory 11. Alternatively, the method may be implemented at the time of acquisition of the CT scan and step 100 may comprises the computer 10 receiving a CT scan transmitted by the CT scanner. The method then comprises a step of building 200 a Finite Element model of the bone, from the obtained CT scan. A Finite Element model (hereinafter, FE model) is a numerical representation of a system or structure, in which the system or structure is discretized into a mesh comprising a plurality of nodes, and three-dimensional elements which are elementary geometrical volumes which corners are defined by the nodes. The elements can be for instance tetrahedron or hexahedron. Furthermore, each element of the FE model may be associated with one or more parameters. In particular, the FE model of the bone may comprise a plurality of elements where each element is associated with mechanical properties derived from the density values of the CT scan. In embodiments, each element of the FE model may be associated with a material type, and mechanical properties may be associated with each material type, such as parameters defining the behavior of the material under load. The building 200 of a FE model comprises a substep 210 of performing segmentation of the CT scan to extract a three-dimensional volume corresponding to the bone. Segmentation is a step of selecting the voxels of the CT scan that correspond to the bone, and excluding the others. The segmentation may be performed by application of a semi-automatic segmentation method, which involves an operator interaction or guidance to perform segmentation. A semi-automatic segmentation method may typically comprise a first rough estimation of an area to be segmented, performed either manually or automatically by application of a first algorithm. Then this rough estimation may be refined during steps involving the operator, such as steps of correction, selection of additional regions of interest, and launching additional image processing algorithms including: - Thresholding, - Contour detection, - Morphologic filters of erosion / dilation, etc. Alternatively, the segmentation may be performed by application of a fully automatic segmentation method, in which the CT scan is processed by a dedicated algorithm for extraction of the three-dimensional volume without any operator intervention. This latter kind of method may be preferred since it reduces the intrinsic variability of the results due to user involvement. An exemplary automatic segmentation method for CT scans is disclosed in US 2021 / 0004956. The building of the FE model then comprises a substep 220 of generating a three- dimensional meshed model of the bone, comprising a plurality of elements, from the three-dimensional volume output by the segmentation step. In other words, the three- dimensional volume output by the segmentation step is meshed into the plurality of elements. This meshing is a step known to the skilled person and may be performed using off the shelf software products such as, for instance, Altair SimLab marketed by the American company Altair. In embodiments, the three-dimensional volume output by the segmentation step, or the three-dimensional meshed model of the bone, may only comprise a portion of the bone. In particular in the case of the femur, the 3D volume or the 3D meshed model obtained from the latter may comprise an upper part of the femur, including the femoral length and a determined length of the femoral shaft. Preferably, the length of the portion of the femur shaft which is represented in the three-dimensional volume is equal to a determined, constant proportion of the total length of the femur. The three-dimensional meshed model of the bone may be defined in accordance with a plurality of parameters such as: - Geometric shape of the elements, which may be for instance selected between tetrahedric and hexahedric, - Size of the elements, - Order number of the elements, where a first order element has nodes only at the corners of the elements, while a second order element also has one node in the middle of each edge. Further, when the bone is a femur, an additional parameter for building the FE model is the length of the portion of femoral shaft, from the femoral head, that is represented in the three-dimensional volume and hence in the FE model, or the proportion of the total length of the femur represented by this length. Indeed, as can be seen in the schematic illustration of figure 5a, the length of the femur shaft has an impact on the lever arm applied by the loads. The building of the FE model then comprises attributing 230 mechanical properties to the meshed model, from the density values of the CT scan. The mechanical properties attributed to the meshed model may comprise: - Mechanical properties associated to the elements of the model and derived from the density values of the CT scan, - A mechanical behavior law, and - A failure criterion. According to a first example, the finite element model may be a model according to the description of [Sas, 2020], wherein the voxel elements of the model are cubic, first order, with an edge size of 2mm. Further, two types of bone materials are defined, cortical and trabecular, such that each element of the FE model is associated to one among the two types of bones. The mechanical behavior law of the two materials is represented by a curve of stress against strain that is represented in figure 4 and comprises: - A first linear phase defined by a Young modulus (E), - A second perfectly plastic phasis defined by the maximal strain (S) and the maximum plastic elongation (εAB) - A third softnening phase defined by a softening module (Ep), and - A fourth perfectly plastic phase defined by the minimal plastic strain (σmin) Thus, the mechanical parameters introduced above E, S, εAB, Ep and σmin are defined for each bone material type, with values provided in the Table 1.

[0002] These parameters are expressed as a function of the bone ash density, which may be defined by the following equation, according to [Keyak, 2005]: ^^^^= 0.0633 + 0.887^^^^Where ^^^^is the calcium hydroxyapatite (CaHA) density which is itself derived from the K2HPO4density by the following equation: ^^^^^= 1.15^^^^^^^– 0.0073 The K2HPO4 density may be inferred from the density values of the CR scan using a K2HPO4 calibration phantom, as is known to the skilled person. The failure criterion in this model is the maximum observed stress. According to another example, the finite element model is a model according to [Duchemin, 2008], wherein the meshing is tetrahedrical of second order, with a volume per element of 1.6 mm3, and a mechanical behavior law is a linear elastic law with a constant Poisson’s coefficient (ν = 0.4) and a Young’s modulus defined according to the hydroxyapatite density of each element, the Young’s modulus being defined by: ^ = 1^^^^^≤ 0 ^ ^ = 10095^^^^^ ^^^> 0 (7) Where ^^^is the hydroxyapatite (HA) density which is derived from the intensity or gray level values associated to the voxels of the CT scan using a calibration phantom as described above. According to still another example, a Finite Element model of a metastatic vertebra may be computed in accordance with [Stadelmann, 2020]. The method then comprises performing numerical simulation 300 of load application on the FE model of the bone, according to pre-determined boundary conditions. The boundary conditions are defined by a set of parameters which may include at least: - Parameters defining the position of load application on the bone, and - Parameters defining a relative orientation of the bone with respect to the direction of the applied loads. With reference to figure 5a is shown an illustration of a configuration of load application on a femur. The arrows represent the position and direction of load application on the femur. The loads are applied vertically on the femoral head, while the other end of the femur is considered fixed. As can be seen from this figure, the parameters involved in the numerical simulation may include, when the bone is a femur, an orientation of the femur in a front plane, a rotation angle of the femur about the femur shaft, also called diaphysis rotation, denoted θ and a position of application of the load on the femoral head. With reference to figure 5b, is shown an illustration of a configuration of load application on a vertebra. As shown in this figure, the vertebral body has two opposite main surfaces and the load may be applied uniformly on the upper main surface, while the bottom main surface is considered fixed. The parameters involved in the numerical simulation may comprise an orientation of the vertebra in a sagittal plane, an orientation of the vertebra in a frontal plane as well as a position of application of the loads on the vertebra. The numerical simulation of load application on the FE model comprises simulating application of an increasing load at the determined point of application of the load, and computing the resulting deformation and strains in the elements of the FE models until the failure criterion is reached. The load value at which the failure criterion is reached is recorded, and it is called simulated failure load. The method then comprises a step 400 of determining a probability of fracture of the bone, and associated uncertainty, for a determined load value applied to the bone, from the result of the numerical simulation, in particular the simulated failure load value determined from the simulation, and from the uncertainties related to the parameters involved in building the finite element models and / or the boundary conditions of the simulation. The probability of fracture of the bone may be determined from a previously established uncertainty on the predictions provided by the numerical simulations, where said uncertainty results from the cumulative uncertainty related to the successive steps of building the FE model and performing the numerical simulation. In particular, uncertainty may result, regarding the building of the FE model, from: - parameters involved in the segmentation, especially if it is a semi-automatic segmentation method (for instance parameters of steps of erosion and dilation), and - parameters involved in the derivation of mechanical properties from the density values of the CT scan. Uncertainty may also result, regarding the numerical simulation, from: - parameters mentioned earlier defining the position of load application, and relative orientation of the bone with respect to the direction of the applied load, and - when the bone is a femur, a length of the portion of femur shaft represented in the FE model. Considering an in-vivo bone subject to a load, let be denoted F the actual rupture load of the bone, ^^the simulated failure load, i.e. the failure load determined from the numerical simulation performed according the method disclosed above, and Δ^^the uncertainty related to the simulated failure load: ^= ^^± Δ^^Evaluating the uncertainty of the failure load value provided by the numerical simulation may be performed by a comparison with the experimental ex-vivo failure load. Thus, the method may comprise a preliminary experimental step 90, wherein, for a plurality of ex-vivo metastatic bones, a comparative experiment is performed by: - performing a CT-scan on an ex-vivo bone, building a CT scan of the bone according to step 200 disclosed above, and performing numerical simulation of load application on the obtained FE model of the ex-vivo bone until the determined failure criterion is fulfilled for a given simulated failure load value, and - performing load application on the ex-vivo bone, according to the same boundary conditions used in the numerical application, to determine an actual failure load. By “same boundary conditions”, it is meant that the position and direction of application of the loads on the bone, as well as the orientation of the bone with respect to said direction, are the same in the numerical simulation and in the experimental load application performed on the ex-vivo bone. It is therefore obtained, for each ex-vivo bone, an experimental failure load value Fexp and a simulated failure load value ^^. From the values obtained for all ex-vivo bones, a standard deviation SD(^^− ^^^^) can be obtained. If one considers that the experimental load application is devoid of uncertainty, i.e. ΔFexp ≈ 0, then: ^ ≈ ^^^^ and ^^^^ ≈ ^^± Δ^^It is further estimated that the average of the differences between the experimental values ^^^^ and the predicted values is null. The uncertainty on the predicted value of the numerical simulation can thus be defined as two times the standard deviation of the differences between experimental values and predicted values, such that: In order to determine relevant uncertainty value for the probability of fracture determined for the in-vivo metastatic bone, all the steps of building of the FE model and performing numerical simulation are identical in the preliminary step 90 and in the steps 200 and 300 performed for the patient. In particular, parameters involved in the segmentation, in the definition of the 3D meshed model, in the derivation of mechanical properties of the FE model, and in the boundary conditions of the numerical simulation, are the same during the comparative experiment 90 and the steps performed for the in vivo bone of the patient. Thus, the uncertainty Δ^^computed at step 90 may be stored in a memory along with the values of all the relevant parameters used during step 90, so that they can be transmitted to a practitioner, who can then use these parameters for performing the same simulation on its patient, and use the uncertainty value to infer a fracture probability from the simulated failure load value he has obtained from the simulation on the patient. A density of probability p(x) of fracture of the bone may then be determined from the computed uncertainty Δ^^as: Where P is the weight of the patient in N, Δ^^is the uncertainty determined at step 90 and ^^is the simulated failure load value determined at step 300. The probability of fracture at a determined load value applied to the bone may then be readily computed from the density of probability, as follows: ^^= 0.5 ∗ ^1 + erf Where VC is the determined load value expressed as a multiple of the patient’s body weight and erf is Gauss error function. In [Bergman, 2001], the force exerted on femur has been measured for a plurality of daily activities, and it has been shown that the maximum force was 3.5 times the body’s weight while descending stairs. From this work, it is thus possible to compute a probability of fracture of a femur in daily activity by computing the probability above for VC = 3.5. On figure 6a is shown an example of density of probability of an actual failure load value, with the probability of fracture at VC=3.5 times the body weight. The colored area represents an example of density of probability of the failure load value simulated using a metastatic femur FE model for a body weight of 1000 N which corresponds to approximately 100 kg, with a simulated failure load value of 10000 N, and an uncertainty of 3 times the body weight i.e.3000 N. The dotted line represents the maximum load value applied on a femur during stairs descent, i.e.35000 N. On figure 6b is shown the probability of fracture of a femur according to the critical value VC, where the curve is the cumulative distribution function, i.e. the probability that the actual failure load value is lower than the critical value VC, with the same example as in figure 6a. This figure 6b shows that the probability is 1.5% when VC equals 3.5 times the patient’s body weight, meaning that the probability that the actual failure load value is lower than 3.5 times the patient’s body weight is 1.5%. Said differently, the femur of the patient has a probability of fracture of 1.5 % during daily activities such as descending stairs. Moreover, the same probability may be computed for various types of activities, such as slow walking, fast walking, standing up, sitting down, etc. In embodiments, the method may thus comprise computing a probability value for a plurality of determined load values corresponding to respective types of activities. The method may then comprise a step 500 of displaying a summary for the practitioner’s review, an example of which is provided in figure 7, in which each type of activity is associated to the probability of risk of fracture of the bone. The activities that are listed in the summary may, in embodiments, be personalized according to the patient’s usual routine. In embodiments, when the metastatic bone is a femur, and the patient has one metastatic femur and a healthy femur, the method may comprise acquiring a CT scan of both femurs and performing the steps 200 and 300 of building a FE model and running numerical simulation to determine a simulated failure load value for both the metastatic femur and the healthy femur. This enables computing a level of weakening of the metastatic femur as compared to the healthy one. The sensitivity of two types of FE models of femur with respect to parameters of segmentation, boundary conditions of the simulation, and mechanical parameters involved in the building of the FE model, was assessed as follows. Materials and Methods Two datasets were used. The Leuven dataset, i.e. the data set of the original model disclosed in [Sas02,2022] consisted of the CT scans (voxel size: 0.4x0.4x0.2 mm3) of 8 femurs with surgical defects, cut at 25 cm (from the top of the head), experimentally tested until failure in a setup mimicking single leg stance. A full description of the experimental setup has been published [Sas02,2020]. A second dataset, so-called Lyon dataset, was used to test replicability of the original model. This dataset consisted of the CT scans (voxel size: 0.7x0.7x1.2 mm3) of 16 intact femurs and 6 femurs (3 pairs) with surgical defects cut at 15 cm (from the lesser trochanter), together with the description of the experimental setup and experimental failure load used. The evaluated workflow (the original model) is a nonlinear finite element model described and validated in [Sas,2020]. The model is originally described as follows. The femur is segmented from the CT scan (voxel size: 0.94x0.94x3 mm3) (no guideline is given on the segmentation technic). Then, the CT volume is resliced in order to get 2 mm isotropic voxels. A first-order hexahedral (brick, voxel-based) mesh is created based on the segmentation, and mechanical properties are assigned from the Hounsfield unit of each voxel to the associated element. The Hounsfield units are first converted to calcium hydroxyapatite (CaHA) density using a phantom of CaHA scanned along with the femur. Any negative density voxel was set to 0 g.cm-3. Surface elements with a CaHA density smaller than 0.25 g / cm3 were removed to account for partial volume effects of the CT scan. Then, the densities were smoothed with a gaussian (SD = 1.27 with limited kernel of 3x3 elements). The CaHA density was then converted into ash density with the relations provided by [Keyak, 2005]. ^^^^= 0.887 × ^^^^^+ 0.0633 The ash density was then converted into mechanical properties according to Table 1 above. In order to accurately reproduce the experimental setup, boundary conditions were carefully applied to the femur model in this study. Specifically, the load was applied on the top surface of a load cup that was placed on the femoral head, while a null displacement was applied to a support embedding the distal diaphysis. These boundary conditions were chosen to closely mimic the experimental conditions, ensuring that the results obtained from the model were as accurate as possible. To solve the model, ParOsol was used, which is an implicit finite element analysis software. This software can simulate linear models. Additionally, a custom modification was applied to the software, as explained in [Flaig, 2011], to enable the simulation of nonlinear models. Elemental results were computed using gaussian quadrature, which is a reduced integration method commonly used in finite element analysis. The model performance was evaluated by comparing the predicted failure load with the experimental ex vivo failure load. For easier comparison, we will use the coefficient of determination (r2), the average error between predicted and experimental failure loads (µ), the standard deviation of errors between predicted and experimental failure loads (SD), and the number of samples (N) to evaluate the model performance. The model performance in the original paper was (Norig._intact = 10, rorig._intact2= 0.90, µorig._intact = -8%, SDorig._intact = 13%) for intact femurs and (Norig._lesion = 10, rorig._lesion2= 0.93, µorig._lesion = 21%, SDorig._lesion = 22%) for femurs with surgical lesions. Some discrepancies between the original and reproduced models are to be noted: - Discrepancies in boundary conditions: the supports were not modelled, the load was directly applied to the top of the femoral head, and the nodes on the surface of the distal diaphysis were applied as a null displacement. While this may have introduced some variability in the results, it has been previously verified that these boundary conditions did not induce any variation on the failure load by modelling the two boundary conditions on three femurs and observing the same failure load for both boundary conditions. - Discrepancies in segmentation technics (no guideline given in the original paper) - Discrepancies in the orientation method of the femurs (Differences between the original and replicated models’ orientation method led to less than 2° of difference in the femur orientation) - Discrepancies in the software used. Whereas the original model used Mimics (Materialise NV) and ParOSol

[0026] , the reproduced model used 3DSlicer (slicer.org) for the segmentation, QCTMA (pypi.org / project / qctma) for the mechanical properties attribution, and Ansys Mechanical (ANSYS, Inc.) for the simulation. These discrepancies were not corrected, as the segmentation and orientation processes were not fully disclosed in the original paper. One major discrepancy was corrected: the original model’s dataset used a CaHA phantom, while the reproduced model used a K2HPO4 phantom, and a conversion law was used to directly get Ash density from K2HPO4 density, while the original model used a conversion law to get CaHA density from K2HPO4 density and then continued with the originally described process. The reproduced model was then corrected upon choosing the conversion method from Leuven’s team. Each sample of Leuven’s and Lyon’s datasets was mechanically tested in a compression test using an Instron 3360 (Norwood, MA, USA) for the Leuven’s dataset, and an Instron 8802 (High Wycombe, England) for Lyon’s dataset. The sample was oriented in a single-leg stance configuration. The sample was preloaded, then small load cycles were applied to remove the nonlinear effects of the initial loading and finally loaded in a quasi-static movement until failure. The distal diaphysis was embedded in polymethylmethacrylate (PMMA) and fixed to the testing machine. A PMMA support was molded to each femoral head and placed between the loading plate and the head to distribute the applied load. For the Lyon’s dataset, a laser scan was used on each sample to get the exact orientation during the experimental tests and minimize the uncertainty of orientation (the surface mesh of the experimental setup was matched to the surface mesh of the segmentation by a point-set registration algorithm). The reproduced model was applied to each sample of the Leuven dataset, and predictions of the original and reproduced models were compared. To evaluate the replicability of the model, the Lyon’s dataset was used. The performance of the reproduced model applied to the Leuven’s dataset was then compared to its performance when applied to the Lyon dataset. To evaluate the performance of models, the mean and standard deviation of the differences between predicted and experimental failure load were used, as were linear regression parameters between predicted and experimental failure loads. The similarity between the original and reproduced model was assessed by the linear regression parameters between the predictions of both models. The value of the differences will be given as a percentage of the relative difference compared to the experimental values. Two global sensitivity analysis were done to evaluate the reproduced model. The objective of the first sensitivity analysis was to evaluate the sensitivity of the predicted failure load with regards to the direct input parameters of the model. The mechanical parameters used to define the non-linear constitutional law are obtained with functions dependent of the mineral density of the bone. In that regard, the objective of the second sensitivity analysis was to evaluate the sensitivity of the predicted failure load with regards to the uncertainties of the conversion laws used to convert the mineral density to the mechanical parameters. The first global sensitivity of the reproduced model was evaluated using the Morris method [Saltelli, 2004] with 10 trajectories and 4 levels for each of the following parameters: femur length, orientation, segmentation, calibration coefficient and diaphysis rotation (Table 2). The femur used for the global sensitivity analysis was chosen among the Leuven dataset. Ranges of variation for each parameter were chosen by repeatedly applying the nominal values and then computing the standard variation of each measured parameters. The range corresponds to the nominal value + / - two standard deviations. Thus, the results of the Morris method account mostly for the uncertainty due to the operator. The elementary contribution of each parameter on the variation of the failure load will be presented as a percentage of the failure load obtained with the nominal values. Table 2. Parameters studied with the Morris method. Nominal Parameters value Min Max Parameter description Femur length (mm) 250 230 250 Working length of the femur Orientation (°) 12 6 18 Orientation of the diaphysis axis in the frontal plane 0 : reference segmentation -1 : erosion with 1 voxel kernel Segmentation 0 -1 +2 from the reference +1 an +2 : dilation with respectively 1 and 2 voxels kernels from the reference Calibration coefficient 8.57 ×-4 -4 -4Slope of the Hounsfield units to (g.cm-3) 10 8.47 × 10 8.67 × 10 K2HPO4 conversion equation Diaphysis rotation (°) 0 -8 8 Rotation around the diaphyseal axis The second global sensitivity of the reproduced model was evaluated using the Morris method with 5 trajectories and 4 levels for each parameter described in Table 3, totaling 70 simulations. The femur used for the global sensitivity analysis was the same as the one used for the first sensitivity analysis. Ranges of variation for each parameter were chosen by computing the uncertainty of the coefficient of regressions with respect to the method described in [Taylor, 1997]. The uncertainty corresponds to one standard deviation. The range for each parameter used in the sensitivity analysis corresponds to the nominal value + / - the uncertainty. The sensitivity analyses were conducted by using the SALib Python package (https: / / github.com / salib / salib). Table 3 Conversion laws of mechanical parameters. Their nominal value and uncertainty were computed from raw data extracted from the article where they were taken from. Conversion Equation Parameter Nominal Incertainty law value c_Ee 14900 320 ϱAsh-> E E = c_Ee * ϱAsh^ d_Ee d_Ee 1.86 0.135 c_S 102 4.4 ϱAsh-> S S = c_S * ϱAsh^ d_S d_S 1.8 0.201 εAB= c_epsAB_cort * ϱAsh^ c_epsAB_cort -0.01 0.02998 d_epsAB_cort d_epsAB_cort 0.0184 0.003544 ϱAsh -> εAB εAB= c_epsAB_trab * ϱAshc_epsAB_trab 0.0241 0.02998 ^ d_epsAB_trab d_epsAB_trab 0.00189 0.003544 Ep= Epp_cort Epp_cort -1000 1304 ϱAsh-> EpEp= c_Epp_trab * ϱAsh^ c_Epp_trab -2080 1304 d_Epp_trab d_Epp_trab 1.45 0.34 σ = c_sigma_ c_sigma_min 43.1 41.43 ϱ min min * ϱAsh Ash -> σmin ^ d_sigma_min d_sigma_min 1.81 0.54 ϱK2HPO4 : K2HPO4 density ϱCHA : CHA density ϱAsh : Ash density E : Young’s modulus S : Ultimate stress εAB : Max plastic elongation Ep : Softening module σmin : Min stress Results There is a high correlation between the original and reproduced models (rreproduction2= 0.95). Dispersions of errors are similar (SDoriginal = 15%, SDreproduced = 16%); however, the systematic bias is doubled for the reproduced model (µoriginal = 11%, µreproduced = 25%). Correlations with experimental data are similar (roriginal2= 0.96, rreproduced2= 0.98). Therefore, the model is reproducible, despite having a 14% additional systematic error. The performance of the reproduced model shown previously using the Leuven’s dataset tends to degrade when using the Lyon’s dataset (rLyon2= 0.60, µLyon = -39%, SDLyon= 28%). These results contain the performance for intact femurs only (rintact2= 0.66, µintact= -50%, SDintact= 17%) and femurs with surgical lesions (rlesion2= 0.03, µlesion= -11%, SDlesion= 32%). There is little to no correlation at all with the experiment for this group. The systematic bias is smaller for femurs with lesions and larger for intact femurs. The dispersion of errors, however, is better for intact femurs and worse for femurs with lesions. Five simulations (spread on 2 trajectories) out of 60 crashed and could not give any result. The 3 defective trajectories were removed from the calculation, leading to the final Morris analysis having 8 trajectories relying on 48 simulations. All the crashed simulations were simulations with the maximum values for the calibration coefficient and the segmentation. The results are in percentage of the failure load obtained with the nominal values (F = 3766 N). The Morris method shows that the calibration coefficient is the most influent parameter (mean influence of 12.5 % on the failure load). The orientation has a mean influence of 10.1 %, the segmentation has a mean influence of 8.6 %, the length has a mean influence of 6.6 % and the diaphysis rotation has a mean influence of 4.4 %. The Morris analysis of the parameters of the conversion laws shows influences superior to 5% on the failure load of 6 parameters out of 13. The parameters of the ultimate stress (d_S and c_S) have a mean influence of respectively 44 % and 14 %. The parameters of the minimum stress (d_sigma_min and c_sigma_min) have a mean influence of respectively 26 % and 7 %. The parameter of the softening modulus of the cortical (Epp_cort) has a mean influence of 12 % and finally, one parameter of the plastic max elongation (c_epsAB_cort) has a mean influence of 10 %. REFERENCES [Eggermont, 2022] Eggermont and al.: “A Patient-Specific Fracture Risk Assessment Tool for Femoral Bone Metastases: Using the Bone Strength (BOS) Score in Clinical Practice”, Cancers, 2022, 14, 5904 [Sas, 2022] Sas A, Tanck E, Wafa H, van der Linden Y, Sermon A, van Lenthe GH, “Fracture risk assessment and evaluation of femoroplasty in metastatic proximal femurs. An in vivo CT-based finite element study.” J Orthop Res. 2022; 1-10. doi:10.1002 / jor.25331 [Sas02,2022] Sas, A.; Sermon, A.; van Lenthe, G.H. Experimental Validation of a Voxel-Based Finite Element Model Simulating Femoroplasty of Lytic Lesions in the Proximal Femur. Sci Rep 2022, 12, 7602, doi:10.1038 / s41598-022-11667-x. [Sas, 2020] Sas, A. et al., “Nonlinear voxel-based finite element model for strength assessment of healthy and metastatic proximal femurs”. Bone reports, 12? 100263, 2020, https: / / doi.org / 10.1016 / J.bonr.2020.100263 [Sas02, 2020] Sas, A.; Van Camp, D.; Lauwers, B.; Sermon, A.; van Lenthe, G.H. Cement Augmentation of Metastatic Lesions in the Proximal Femur Can Improve Bone Strength. Journal of the Mechanical Behavior of Biomedical Materials 2020, 104, 103648, doi:10.1016 / j.jmbbm.2020.103648. [Keyak, 2005] Keyak, J.H.; Kaneko, T.S.; Tehranzadeh, J.; Skinner, H.B. Predicting Proximal Femoral Strength Using Structural Engineering Models. Clin Orthop Relat Res 2005, 219–228, doi:10.1097 / 01.blo.0000164400.37905.22. [Duchemin, 2008] Duchemin, L.; Mitton, D.; Jolivet, E.; Bousson, V.; Laredo, J.D.; Skalli, W. An Anatomical Subject-Specific FE-Model for Hip Fracture Load Prediction. Computer Methods in Biomechanics and Biomedical Engineering 2008, 11, 105–111, doi:10.1080 / 10255840701535965 [Stadelmann, 2020] Marc A. Stadelmann, Denis E. Schenk, Ghislain Maquer, Christopher Lenherr, Florian M. Buck, Dieter D. Bosshardt, Sven Hoppe, Nicolas Theumann, Ron N. Alkalay, Philippe K. Zysset, “Conventional finite element models estimate the strength of metastatic human vertebrae despite alterations of the bone's tissue and structure”, PII: S8756-3282(20)30378-1 DOI: https: / / doi.org / 10.1016 / j.bone.2020.115598 Reference: BON 115598 [Bergman, 2001] Bergman et al. ‘Hip contact forces and gait patterns from routine activities”, in J. Biomech., 2001 [Flaig, 2011] Flaig, C.; Arbenz, P. A Scalable Memory Efficient Multigrid Solver for Micro-Finite Element Analyses Based on CT Images. Parallel Computing 2011, 37, 846–854, doi:10.1016 / j.parco.2011.08.001. [Saltelli, 2004] Saltelli, A. (Andrea) Sensitivity Analysis in Practice : A Guide to Assessing Scientific Models; Chichester ; Hoboken, NJ : Wiley, 2004; ISBN 978-0- 470-87093-8. [Taylor, 1997] An Introduction to Error Analysis: The Study of Uncertainties in Physical Measurements.

Claims

CLAIMS 1. A computer-implementing method for determining a risk of fracture of a metastatic bone of a patient, consecutive to a load being applied to the metastatic bone, wherein the bone is either a femur or a vertebra, the method comprising: - Obtaining (100) a Computed Tomography (CT) scan of the bone, the CT scan comprising a plurality of voxels associated to respective density values, - Building (200), from the CT scan of the bone, a Finite Element Model, comprising: o segmenting (210) the CT scan to extract a three-dimensional volume corresponding to the bone, o generating (220) a three-dimensional meshed model of the bone, comprising a plurality of elements, from the three-dimensional volume, o attributing (230) mechanical properties to the meshed-model, said mechanical properties being derived from the density values of the CT scan, - performing numerical simulation (300) of load application on the FE model of the bone, according to pre-determined boundary conditions, to determine a simulated failure load value, and - determining (400) a probability of fracture of the bone, for a determined load value applied to the bone, based on the simulated failure load value and a determined uncertainty related to the building of the FE model and the numerical simulation.

2. The method according to claim 1, wherein the boundary conditions of the numerical simulation comprise parameters defining the position of load application, and a relative orientation of the bone with respect to the direction of the applied loads.

3. The method according to any of the preceding claims, comprising a preliminary step of determining an uncertainty value related to the building of the FE model and the numerical simulation, comprising: - For a plurality of ex-vivo bones: o Performing a CT-scan on the ex-vivo bone, building a Finite Element model of the ex-vivo bone, and implementing a numerical application of load application on the FE model of the ex-vivo bone until a predetermined failure criterion is fulfilled for a simulated failure load, the numerical simulation of load application being performed according to determined boundary conditions on the position of load application and the relative orientation of the bone with respect to the direction of the applied loads, o Performing load application on the ex-vivo bone according to the same boundary conditions as used in the numerical application, to determine an actual failure load, and o Computing the uncertainty value from the standard deviation between the simulated failure loads and the actual failure loads for the plurality of ex-vivo bones.

4. The method according to any of the preceding claims, wherein determining a probability of fracture of the bone for a determined load value applied to the bone, from the simulated failure load value and the determined uncertainty value, comprises: - Determining, from the determined uncertainty value and from the simulated failure load value obtained from the numerical simulation performed on the bone, a density of probability of the failure load, and - Computing, for the determined load value applied to the bone, a probability of fracture as the integral of the density of probability up to the determined load value.

5. The method according to any of the preceding claims, wherein the determined load value is a multiple of the body weight of the patient.

6. The method according to any of the preceding claims, wherein the segmentation is performed by implementation of an automatic segmentation algorithm.

7. The method according to any of the preceding claims, comprising computing a probability of fracture of the bone, and associated uncertainty, for each of a plurality of determined load values corresponding to respective average load values applied to the bone in a plurality of situations or activities to which the patient may be subject, and displaying, for each situation or activity, the corresponding probability of fracture.

8. The method according to any of the preceding claims, wherein the bone is a metastatic femur, and the method is applied to the metastatic femur of a patient and to the other femur of the patient, and further comprises displaying, for each femur of the patient, the probability of fracture of each femur for a same applied load value.

9. A computer-program product, comprising code instructions for implementing the method according to any of the preceding claims, when it is implemented by a computer.

10. A device for determining a risk of fracture of a metastatic bone of a patient, consecutive to a load being applied to the metastatic bone, wherein the bone is either a femur or a vertebra, the device comprising at least a computer and a memory, and being characterized in that it is configured for implementing the method according to any of claims 1-9.

11. The device according to claim 10, further comprising a CT-scanner in communication with the computer, the CT-scanner being configured to acquire CT-scans and to communicate the acquired CT-scans to the computer.

Images