Computational tibial bone remodeling over a population after total knee arthroplasty: A comparative study

Abstract Periprosthetic bone loss is an important factor in tibial implant failure mechanisms in total knee arthroplasty (TKA). The purpose of this study was to validate computational postoperative bone response using longitudinal clinical DEXA densities. Computational remodeling outcome over a population was obtained by incorporating the strain‐adaptive remodeling theory in finite element (FE) simulations of 26 different tibiae. Physiological loading conditions were applied, and bone mineral density (BMD) in three different regions of interest (ROIs) was considered over a postoperative time of 15 years. BMD outcome was compared directly to previously reported clinical BMD data of a comparable TKA cohort. Similar trends between computational and clinical bone remodeling over time were observed in the two proximal ROIs, with most rapid bone loss taking place in the initial months after TKA and BMD starting to level in the following years. The extent of absolute proximal BMD change was underestimated in the FE population compared with the clinical subject group, which might be the result of significantly higher initial clinical baseline BMD values. Large differences in remodeling response were found in the distal ROI, in which resorption was measured clinically, but a large BMD increase was predicted by the FE models. Multiple computational limitations, related to the FE mesh, loading conditions, and strain‐adaptive algorithm, likely contributed to the extensive local bone formation. Further research incorporating subject‐specific comparisons using follow‐up CT scans and more extensive physiological knee loading is recommended to optimize bone remodeling more distal to the tibial baseplate.

tibiae. Physiological loading conditions were applied, and bone mineral density (BMD) in three different regions of interest (ROIs) was considered over a postoperative time of 15 years. BMD outcome was compared directly to previously reported clinical BMD data of a comparable TKA cohort. Similar trends between computational and clinical bone remodeling over time were observed in the two proximal ROIs, with most rapid bone loss taking place in the initial months after TKA and BMD starting to level in the following years. The extent of absolute proximal BMD change was underestimated in the FE population compared with the clinical subject group, which might be the result of significantly higher initial clinical baseline BMD values.
Large differences in remodeling response were found in the distal ROI, in which resorption was measured clinically, but a large BMD increase was predicted by the FE models. Multiple computational limitations, related to the FE mesh, loading conditions, and strain-adaptive algorithm, likely contributed to the extensive local bone formation. Further research incorporating subject-specific comparisons using followup CT scans and more extensive physiological knee loading is recommended to optimize bone remodeling more distal to the tibial baseplate. interventions, but despite reduced revision rates, the number of primary TKA failures is increasing as a result of the aging population and the acceptance of TKA in younger patients. 1 Two common causes of long-term implant failure are aseptic loosening and periprosthetic fracture, 1,2 which are linked to stress shielding-related bone loss as observed in longitudinal DEXA studies, [3][4][5] in line with Wolff's law.
Clinical DEXA studies typically display a significant spread in bone density changes, attributed to various sources of variation in patient characteristics, such as differences in initial bone mineral density (BMD), preoperative knee alignment, and subject body mass index (BMI).
For instance, initial mediolateral (ML) bone density distributions may depend on native knee alignment, since clinical studies have demonstrated a higher medial baseline BMD in varus knees, while valgus knees typically have greater initial lateral bone density. 3,[6][7][8] Subsequently, knee alignment may be changed during surgery, leading to changes in joint loads and load transfer to the periprosthetic bone. Hence, a mechanical TKA alignment may lead to relatively more medial bone loss in constitutional varus knees, and more laterally concentrated bone loss in native valgus knees.
This trend has been observed in numerous clinical DEXA studies. 3,4,6,[8][9][10] A higher initial BMD has also been associated with greater (relative) proximal bone loss regardless of knee alignment, by a computational and a clinical study, 6,11 respectively, causing mean proximal BMD to converge to a fixed density range after 2 years.
Furthermore, clinical studies have demonstrated a positive correlation between subject BW measures and proximal BMD levels several years after TKA. 8,9,12,13 Conversely, no pronounced effect of age and sex on tibial bone loss has been reported. Since the age range of a primary TKA cohort is typically limited, different TKA studies have been unable to demonstrate age-related BMD decline in the preoperative tibia and subsequent remodeling. 9,14 However, age-related bone loss is generally more pronounced in postmenopausal women, 15 accounting for higher baseline BMD levels found in male TKA patients compared with female patients. 14 No significant differences in postoperative density changes were found by sex in various studies. 8,12,14,16,17 One study reported significantly less bone loss in lateral and distal regions in male patients, 18 potentially due to corrective alignment change related to preoperative varus deformity, as constitutional varus is more common in men than in women. 19 Investigating the effect of various sources of variation on periprosthetic bone changes in a clinical setting in more detail would require long-term follow-up studies with large patient cohorts. An alternative way to gain more understanding about the relative effects of these parameters is through computational modeling. Previous finite element (FE) models have assessed periprosthetic tibial bone loss using strain differences, 20,21 and by subsequent modeled bone loss through strain-adaptive bone remodeling. 11,22 Current strainadaptive remodeling theories have been established and refined based on femoral bone changes following total hip replacement, 23,24 but to our knowledge have not been validated before against clinical outcome in the tibia. In the current study, computational bone remodeling outcome in a TKA cohort was compared against results of a longitudinal clinical DEXA study in a different patient group. [25][26][27] The computational results were furthermore used to investigate the relative effects of patient characteristics on periprosthetic bone remodeling in more detail.

| MATERIALS AND METHODS
FE models were created using a custom-made workflow, 11 based on a Japanese lower limb CT data set of 26 tibiae from 14 subjects who were scanned prior to TKA surgery. The model setup consisted of the following consecutive steps: (a) CT scan processing, (b) FE mesh generation, (c) material property assignment, and (d) application of boundary and loading conditions.
In the initial processing step, knee alignment angles were measured and tibiae were segmented from the available CT scans, taken in supine position. The hip-knee-ankle (HKA) and tibial varus-valgus (VV) angles of each knee were measured in the anteroposterior (AP) view using CT scan annotations in Slicer 3D 28 ; the HKA angle was defined as the angle between the femoral and tibial mechanical axes, while the VV angle was the angular offset of the joint line perpendicular to the tibial mechanical axis. Positive knee angles were directed toward varus alignment. Based on the measured knee angles, each knee was allocated to be in varus, neutral or valgus alignment.
Knees were considered to be in varus if the HKA angle was greater than or equal to 3 . Subdivision between neutral and valgus knees was made based on the measured tibial VV angle; knees were assigned to be in valgus in case of a valgus tibial joint line (VV < 0 ).
Subject characteristics including the measured knee angles and resulting alignment distributions are indicated in Table 1.
The tibiae were automatically segmented based on boundary enhancement filtering and graph cut optimization 29 ; the resulting segmentation was manually adjusted using Slicer 3D in case incorrect local bone edges were detected. Surface meshes were generated from the obtained binary voxel masks 30 and smoothed using curvature flow. 31 The bones were subsequently aligned according to the mechanical axis, 32 with the largest inertial axis being defined as the longitudinal axis and neutral internal rotation referencing the medial third of the tibial tubercle. 33 T A B L E 1 Subject characteristics indicated by mean (range) to the posterior condyles, which was adopted to compensate for flexion and extension gaps. 34 A posterior tibial slope of 5 was adopted, based on the general surgical recommendation for cruciate-retaining implants to match the patient's anatomy up to 5 of posterior slope, and an average anatomical slope of 10 and 6 in neutral and varus knees, respectively, in a Japanese population. 35 The tibial resection level was defined 8 mm distally from the lowest point of the highest condyle, and the internal/external rotation of the RP tray was optimized to maximize the coverage of the resected bone surface. The correct implant size was set to be the largest tray size which could be placed on the resection surface with a maximum overhang below 2 mm, in line with a reported tibial coverage study. 36 Bone coverage, defined as the relative resection plateau surface area covered by the base plate, was computed to numerically assess the achieved implant position and indicated that realistic implant positions were achieved, since values were in line with results of previous clinical and computational studies. 13,36,37 Established implant sizes and bone coverage ratios are also indicated in Table 1. In the preoperative models, the COPs of the forces, according to the subject's native deformity, were connected to the closest nodes on the proximal tibial surface using springs; the number of connected nodes was determined as function of the total related contact area, and spring constants were individually assigned based on distance and a compressive modulus of 9 MPa representing the intermediate articular cartilage. 44 Neutral preoperative alignment was considered to be consistent with the planned mechanical implant alignment. In native varus alignment, a 3 HKA angle and a 5 VV angle were adopted in the preoperative situation, based on the constitutional varus HKA alignment found over multiple cohorts, 19,45,46 in combination with an additional 2 tibial varus offset in the anatomical joint line. 47 Valgus knees were represented preoperatively using a neutral 0 HKA angle in combination with a À3 tibial VV angle.
The averaged strain energy density (SED) after application of the three activity peak loads was considered as measure for bone strains during daily living. Subsequent iterative bone density changes were simulated using strain adaptive remodeling, with the difference between local preoperative and postoperative SED per unit bone mass, S ref and S, respectively, considered as stimulus for density change in time dρ=dt ½ . 23 If the relative local difference was lower than 35%, the stimulus fell into a lazy zone and no net remodeling was assumed. Outside of this range, the rate of local bone apposition or resorption was dependent on its available free bone surface a, representing the porosity and specific surface and determined based on the corresponding bone density ρ. 48 Bone associated with greater free surface density a was assumed to be more responsive to changes in SED, since remodeling activity takes place at these free surfaces.
Definition of the local bone remodeling rate dρ=dt ½ following the strain adaptive theory, 10   F I G U R E 2 Schematic superior view of medial and lateral COP positions relative to the tibial tray during activity peak For Peer Review loads; markers are scaled based on the extent of the corresponding forces indicated in Table 2 at a mean of 11.4 years, 27

| RESULTS
A comparison between available preoperative parameters in the different Japanese TKA populations used in computational and clinical remodeling is shown in Table 3. No significant differences in patient demographics were encountered, but preoperative medial ROI BMD was found to be significantly higher in the clinical DEXA scans, taken 2 weeks prior to TKA, than in the FE models. Other FE subject characteristics were tested against initial BMD, and against relative and absolute 15-year ROI BMD difference ( Note: Quantitative subject characteristics (mean ± SD) were assumed to be normally distributed in all subject groups, and were tested against difference using the two-tailed Z-test.

| DISCUSSION
In the current study, the results of tibial bone remodeling simulations were compared against a longitudinal clinical DEXA study in a comparable patient group. [25][26][27] The computational results were also used to   The rate of computational bone remodeling over postoperative time indicated that the vast majority of tibial density changes took place within the first 2-3 years, with the highest remodeling rates found in the initial 6 months after surgery (Figure 4), in line with the bone changes in the clinical population. 25 Other DEXA studies also reported the greatest bone loss to occur between initial follow-up time points within 2 years, but found ongoing density changes after  (Table 4).
Computational remodeling was performed on a Japanese TKA population different from the Japanese patients in the clinical remodeling study, [25][26][27] since no preoperative CT scans of the DEXA study were available. No significant differences were encountered in patient demographics of both populations (Table 3). However, preoperative BMD was significantly lower in the medial ROI over the FE tibiae.  Figure 4 and Table 4.
Differences in preoperative BMD between both populations ( Table 3) and underestimation of the net bone loss in the computational group (Table 4) could be the result as well of previously reported effects of (potential) differences in patient characteristics (Table 3). Since preoperative joint angles of the clinical population were not reported, it could be that native varus deformation was more prevalent over these subjects, leading to relatively more initial densification in the medial condyle, 6,54,55 generally increased density over the preoperative proximal tibia, 6 and increased medial bone loss following TKA. 3,6-8 Considering the average age difference between FE and clinical subjects, it could also be that osteoarthritis (OA) was more progressed over the older clinical TKA population; knee OA has also been related to increased constitutional varus angles and higher local proximal baseline BMD, respectively, 46,56 although age was not found to be directly related to the knee angles measured in the FE population (Table 5). Higher proximal preoperative BMD was related to increased bone loss in a previous clinical study, 6 which was in line with the differences in net bone loss between the clinical and the FE group (Table 4), and with the correlations between baseline BMD and 15-year BMD change within the FE group (Table 5).
Despite reported in numerous clinical studies, 3,6-8 no relation was found between measured knee angles and baseline ML density values in the current FE subjects. This could be due to the fact that the alignment angles were measured in supine, non-weight-bearing CT scans, which may differ from alignment in a weight-bearing position. 57 On the other hand, the difference in sex distribution over the assigned alignment subgroups, with 71% of male knees considered in preoperative varus versus 42% of female knees, was in line with the finding that constitutional varus is more common in men than in women, 19 suggesting a reasonable subdivision of the knees over the alignment groups was made.
The effect of tibial alignment on postoperative remodeling was also not significant in the current study, in contrast to a previous computational study. 11 23,24 in which extensive bone formation due to excessive local strain concentration was not reported. As a result, the current remodeling algorithm was mainly derived from observed femoral bone loss as response to decrease in local strains, and behavior of bone formation was currently assumed to be inversely related to bone resorption. However, net bone resorption was found to occur at a much higher rate than bone (re)formation following long-term changes in mechanical loading, 60,61 suggesting a higher lazy zone threshold and decreased sensitivity to be used in case of local SED increase, and additional interactions could take place. For instance, post-yield bone behavior was not accounted for in the bone response and material properties, while an average compressive strain increase of 6.9 and 11.6% in spongy bone of female proximal tibiae were reported for yield and ultimate failure strains, respectively. 62

| CONCLUSIONS
Based on the comparison between computation bone remodeling outcome and clinical density changes over similar populations, we can conclude that the current strain-adaptive remodeling algorithm is able to predict the course of bone density changes over time in the proximal tibial ROIs, but not in the distal ROI. Extensive distal bone formation is likely caused by simplifications in implant fixation, (distal) loading conditions and strain-adaptive theory. To improve computational remodeling, it is recommended to perform further research using intrasubject comparisons, based on longitudinal clinical CT data and physiological load cases including soft tissue representations.
Being able to reliably predict tibial periprosthetic remodeling is helpful to guide clinical practice, by identifying risk factors in implant design, surgical technique and tibial features for potential long-term failure due to excessive regional bone loss.

ACKNOWLEDGMENTS
This study was funded by DePuy Synthes Joint Reconstruction, Leeds, UK. One of the authors is an employee of DePuy Synthes.

DATA AVAILABILITY STATEMENT
Data available on request due to privacy/ethical restrictions.