Multiple Subchondral Bone Cysts Cause Deterioration of Articular Cartilage in Medial OA of Knee: A 3D Simulation Study

Aims To investigate the impact of subchondral bone cysts (SBCs) in stress-induced osseous and articular variations in cystic and non-cystic knee models using finite element analysis. Materials and Methods 3D knee joint models were reconstructed from computed tomography (CT) and magnetic resonance imaging (MRI). Duplicate 3D models were also created with a 3D sphere mimicking SBCs in medial tibia. Models were divided into three groups. In group A, a non-cystic knee model was used, whereas in groups B and C, SBCs of 4 and 12 mm size were simulated, respectively. Cyst groups were further divided into three sub-groups. Each of sub-group 1 was composed of a solitary SBC in the anterior half of tibia adjacent to joint line. In sub-group 2, a solitary cyst was modeled at a lower-joint location, and in sub-group 3, two SBCs were used. All models were vertically loaded with weights representing double- and single-leg stances. Results During single-leg stance, increase in subchondral bone stress in sub-groups B-1 and B-3 were significant (p = 0.044, p = 0.026). However, in sub-group B-2, a slight increase was observed than non-cystic knee model (9.93 ± 1.94 vs. 9.35 ± 1.85; p = 0.254). All the sub-groups in group C showed significantly increased articular stress (p < 0.001). Conversely, a prominent increase in peri-cystic cancellous bone stress was produced by SBCs in groups B and C (p < 0.001). Mean cartilage shear stress in sub-groups B-1 and B-2 (0.66 ± 0.56, 0.58 ± 0.54) was non-significant (p = 0.374, p = 0.590) as compared to non-cystic model (0.47 ± 0.67). But paired cysts of the same size (B-3) produced a mean stress of 0.98 ± 0.49 in affected cartilage (p = 0.011). Models containing 12 mm SBCs experienced a significant increase in cartilage stress (p = 0.001, p = 0.006, p < 0.001) in sub-groups C-1, C-2, and C-3 (1.25 ± 0.69, 1.01 ± 0.54, and 1.26 ± 0.59), respectively. Conclusion The presence of large-sized SBCs produced an increased focal stress effect in articular cartilage. Multiple cysts further deteriorate the condition by increased osseous stress effect and high tendency of peripheral cyst expansion in simulated cystic knee models than non-cystic knee models.


INTRODUCTION
Osteoarthritis (OA) of the knee joint is a common degenerative joint disease. The advanced cases of OA are usually associated with formation of subchondral bone cysts (SBCs) (Ondrouch, 1963). Nearly 50% of knee OA cases present with SBCs (Wu et al., 2007;Raynauld et al., 2008). In the majority of cases, these are spherical or ellipsoidal in shape and located at the articular surface. These are thought to be specific bony adaptations in the weight-bearing areas of joints (Crema et al., 2010;Javaid et al., 2010). Two different theories (synovial intrusion theory and bony contusion theory) described the formation of SBCs in weight-bearing regions (Freund, 1940;Landells, 1953;Rhaney and Lamb, 1955;Ondrouch, 1963). However, more recently the concept of vascular pathology (hypertension) has been proposed about the formation of SBCs in non-load bearing regions (Chan et al., 2017).
Though there is no consensus about the exact relationship between SBCs and knee pain, it is widely accepted and several recent studies have argued that there is an association between SBCs and knee pain in OA knee patients (Landells, 1953;Kornaat et al., 2006;Torres et al., 2006;Javaid et al., 2012). In clinical practice, X-rays remain the key investigation for the diagnosis of osteoarthritis. Advanced imaging modalities, such as CT and MRI, show more accurate and precise structural morphology of these cystic lesions Roemer et al., 2009;Crema et al., 2010). In a clinical study using MRI, Raynauld et al. (2008) found the correlation between change in mean cyst size in millimeters and loss of medial femoral condyle cartilage in individuals with knee osteoarthritis. The exact relationship between bone marrow lesions (BMLs) and SBCs is not clear; however, some recent studies have advocated that BMLs may develop into SBCs (Carrino et al., 2006;Crema et al., 2008Crema et al., , 2010. On the other hand, biomechanical studies addressed stress changes in the osseous elements (Durr et al., 2004;McErlain et al., 2011). A 2D finite element (FE) study demonstrated that the SBCs produce stress effect and stressinduced micro-fractures in the femoral head. The authors hypothesized that these micro-fractures might be the first step in the development of SBCs in OA (Durr et al., 2004). Another study advocated that SBCs produce increased osseous stress around the cyst (McErlain et al., 2011). To our best knowledge, all of these studies have focused on stress changes in osseous elements only while ignoring subchondral bone and articular cartilage degeneration (Durr et al., 2004;McErlain et al., 2011).
The detailed effect of multiple cysts on the knee joint is lacking in the literature. Moreover, a comparison between a single cyst and multiple cysts in relation to cyst location has not been described before. It is meaningful to investigate associations between subchondral bone cyst parameters (e.g., number, size, and location) and associated articular degeneration in OA pathology. Finite element analysis is well known and widely accepted in orthopedics (Bosiakov et al., 2017;Trad et al., 2018;Bini et al., 2019). In this technique, 3D models are used with simulated loads to analyze and predict the outcomes. The aims of this simulation study were to (a) establish the biomechanical stress-strain profile of different-sized SBCs in tibia; (b) evaluate any correlations between cyst number, size, position, and quantitative stress changes; and (c) determine intra-osseous, intra-cystic compression variations and the resultant stress changes in subchondral bone and articular cartilage using the non-cystic and cystic 3D knee models. The results of this simulation study provide insight into the role of SBCs on biomechanical aspect of articular degeneration in OA pathogenesis.

Three-Dimensional Modeling
Institutional ethical approval was obtained before this study, and it was performed according to the Declaration of Helsinki. CT and MR images of an individual (52 years old, male) with early medial compartment knee OA were used for modeling. There was no past history of trauma. The geometrically accurate three-dimensional (3D) bony model of the knee joint was generated from CT scan images of the right knee in the neutral unloaded position. The 1.5 mm slice thickness CT DICOM images were imported into Mimics 10.1 software (Materialise, Leuven, Belgium) to reconstruct the surface geometry of the femur, tibia, and fibula. MR images were used for the reconstruction and reference points of soft tissues including the ligaments (ACL, PCL, MCL, and LCL), menisci and cartilaginous components of the distal femur, tibial plateau, and tibio-fibular joint (TFJ). In Geomagics Studio 11.0 software (Raindrop Company, United States) IGS files of the bony elements were used to get individual structure's volumetric form.     Soft tissues were also handled and assembled in Geomagics. These models were then assembled and exported as STP files representing the 3D knee complex as shown in Figure 1A.
The coordinate axes of the assembled models were aligned as the X axis pointed medially (lateral to medial femoral condyle), the Y axis pointed posteriorly (anterior to posterior), and the Z axis pointed upward (knee to hip). To simulate the SBCs of different size, number, and location, computeraided design (CAD) software (ProE CREO 3.0 PTC Corp., United States) was used. The 3D model of SBC is depicted in Figure 1B. To illustrate the exact fit SBC-bone model, there was no gap between bone and cyst capsule. The bone cysts were modeled with spherical morphology. Two different-sized (diameters: 4 and 12 mm), single and paired cysts located in upper and slightly lower position in relation to joint line were simulated (see Figure 2 and Table 2 for cysts parameters). In total, seven knee models were simulated in this study. All the models were meshed using the software HyperMesh 14.0 (Altair Engineering, Inc., United States). Based on mesh convergence study, the mesh size of 0.3 mm for bones, soft tissues (ligaments, menisci, cartilages) and SBCs was used. The type of mesh used in this study was quadratic tetrahedral elements to minimize elements stiffness. The cortical, cancellous, and subchondral components of bones were established as shown in Figure 1A.

Finite Element Modeling and Material Properties
The reconstructed 3D models (ING files) were imported into finite element analysis software Abaqus 6.14 (Simulia Corp., United States). Interactions were carried out between different parts of the models using this software. The cortical, subchondral, and cancellous bony portions were assigned the Young's modulus (E) of 7,300, 3,000, and 1,100 MPa and Poisson ratio (υ) of 0.3, 0.3, and 0.26, respectively (Choi et al., 1990;Kim et al., 2011;Qiu et al., 2011;Raja Izaham et al., 2012). FE models of soft tissues included cartilage, menisci, and major ligaments (cruciate and collateral ligaments). The menisci and cartilage were modeled as linear-elastic materials (Donahue et al., 2002). Articular cartilage was defined as Young's modulus of 12 MPa and Poisson's ratio of 0.45, whereas menisci were modeled with Young's modulus of 80 MPa and Poisson's ratio of 0.3 (Hopkins et al., 2010). Interfaces between the cartilage and bones were fully bonded. Contact was assigned between the femoral cartilage and meniscus, between meniscus and tibial cartilage, and between  Frontiers in Bioengineering and Biotechnology | www.frontiersin.org femoral and tibial cartilage on both the medial and lateral sides. Surface friction coefficient of 0.02 was used, which is in the normal range for human articular joints (Mow et al., 1993). The ligament models were defined as rubber-like materials with stress-strain relationship according to previous studies (Takeda et al., 1994;Mesfar and Shirazi-Adl, 2005;Innocenti et al., 2009;Innocenti et al., 2014;Koh et al., 2019).

Validation of 3D Model and Boundary Conditions
A review of the literature shows that in the majority of SBCs cases, a secondary closure is more commonly found than an opening (Landells, 1953). Therefore, we simulated a complete cyst model without any opening at the joint surface. The 3D osseous and soft tissues models and material properties applied in this experiment were assigned according to the previous published studies (see Table 1). The simulated force was applied proximally to the femur, whereas the distal end of the tibia was fixed in all degrees of freedom (McErlain et al., 2011). All the models were tested for two different physiological conditions, i.e., two-leg and one-leg standing. The individual's body weight was 75 kg in this study; therefore, about 50% (375 N) of the body weight was applied to simulate two-leg standing condition and a load of 750 N for single-leg standing in axial direction, as given in Figure 1C.

Grouping
For the ease of results interpretation, models were categorized into three groups. Group A consisted of non-cystic knee model. Group B and group C were composed of SBCs knee models of 4 and 12 mm (diameters), respectively (data from clinical series given in Supplementary Table 1; the smallest SBC was 4 mm and largest 12 mm). We further simulated single and multiple cysts, of different locations (upper and lower) in respect to articular surface, as well as anterior and posterior halves of the medial tibial plateau. Groups B and C were further subdivided into three sub-groups: 1, 2, and 3. In sub-group 1, one cyst was present in the anterior half of medial tibial plateau adjacent to the joint line (Figure 2, represented by B-1 and C-1). In sub-group 2, a solitary cyst was modeled in the lower location of anterior half of medial tibia (∼4 mm from joint line) (Figure 2, represented by B-2 and C-2). Sub-group 3 demonstrated two cysts, one in anterior and another in posterior tibial half (Figure 2, represented by B-3 and C-3). Frontiers in Bioengineering and Biotechnology | www.frontiersin.org

Statistical Analysis
The statistical analysis was done using SPSS 16.0 (SPSS Inc., Chicago IL). Descriptive statistics was used to determine means and standard deviations. One-way ANOVA and multiple comparison least significant difference (LSD) tests were used to determine the mean difference. p-value was considered significant when p < 0.05. Nodal analysis was done in the region of interest (ROI), using selected nodes (N = 20) from each sub-group of 3D model. Von Mises stress (VMS), Tresca stress, and cartilage shear stress distributions in the subchondral bones, the cancellous bones, and cartilage were analyzed. Strain analysis of the SBCs was also calculated.

RESULTS
The total number of nodes and elements in the non-cystic knee model were 1,034,891 and 4,911,748 respectively, whereas the total number of nodes in solitary small cyst

Von Mises Stress Distribution in Subchondral Bones
The mean VMS in subchondral bone of SBC models was higher than the non-cystic knee model (Figure 3). During single-leg stand condition, the non-cystic model showed the mean stress value of 9.35 ± 1.94 MPa, whereas in group B (smaller-sized SBC), the mean stress values were 37. 71 ± 19.52, 9.93 ± 1.94, and 41.97 ± 20.77 MPa in sub-groups 1, 2, and 3, respectively. There was a significant increase (p < 0.001) in subchondral bone plate stress in all the sub-groups of group C. Further detail is given in Table 3.

Stress Analysis of Peri-Cystic Cancellous Bones
Tresca stress analysis of the cancellous portion of medial tibia showed the least mean stress of 5.27 ± 0.14 MPa (5.01-5.64) in the anterior medial tibial plateau of non-cystic knee model. The bony defect caused by SBCs in group B: 8.56 ± 0.27, 6.52 ± 0.28, and 8.76 ± 0.14 MPa and group C: Frontiers in Bioengineering and Biotechnology | www.frontiersin.org   22.74 ± 1.13, 6.92 ± 0.68, and 23.08 ± 0.86 MPa produced the statistically significant (p < 0.001) stress changes in the individual sub-groups, respectively. The detailed description of stress changes in all the three groups are depicted in Figure 4 and Table 3.

Strain (E) Analysis of Cysts
On applying body weight, there resulted changes in the cysts. The resultant strain in all cysts is summarized in Figure 5 and Table 3. The inner section of individual cysts showed the peak areas of strain. The 12 mm-sized cyst in group C had highest mean strain value of 0.03 in the upper location near the joint line (sub-groups 1 and 3). Same-sized cyst at lower position showed nearly a 5 times less mean strain of 0.006. As the size of cysts increased from 4 to 12 mm, there was a prominent increase in strain (p < 0.001).

Cartilage Shear Stress
The mean shear stress analysis of medial tibial cartilage showed that the presence of a single small-sized cyst (sub-groups 1 and 2 of group B) with mean shear stress of 0.66 ± 0.56 and 0.58 ± 0.54 MPa was not significant as compared to non-cystic model (group A) (p = 0.374 and p = 0.590). However, the paired cysts of same size caused significant increase in stress changes (p = 0.011). On the other hand, the comparison between group A and group C (larger-sized SBCs) showed statistically significant increase in cartilage stress in all sub-groups: C-1, C-2, and C-3 (p = 0.001, p = 0.006, and p < 0.001) with mean stress values of 1.25 ± 0.69, 1.01 ± 0.54, and 1.26 ± 0.59 MPa than noncystic knee model (0.47 ± 0.67 MPa), respectively. These focal increments in medial tibial cartilage are shown in Figure 6a and Table 4. Figure 6b depicts the increase in articular stress during double-leg and single-leg stands.

DISCUSSION
According to our findings, on applying the tested loads (singleand double-leg stance) there was prominent increase in pericystic bone stress. We also noted significant increase (nearly 2 to 3 times) in mean cartilage stress in SBC models as compared to the non-cystic knee model. The mean shear stress of medial tibial cartilage in non-cystic model showed nearly uniform stress distribution (represented by letter "A" in Figure 6a), and the peak stress was in the antero-medial part of articular cartilage because in the standing position the medial compartment of the knee bears up to 2.2 times as much loading as the contrary lateral compartment (Johnson et al., 1980;D'Lima et al., 2007;Kutzner et al., 2010). However, SBCs increased the cancellous and subchondral bone stress, and there was a focal increase in the articular shear stress. For the cysts in group B, sub-group 3 affected the larger area in both the osseous and cartilaginous components. These altered areas of stress caused significant increase in cartilage stress. The small-sized solitary cyst (subgroups 1 and 2 of group B) had higher mean cartilage stress than the non-cystic model but were statistically non-significant.
On the other hand, there was a significant increase in focal cartilage stress (which may cause degeneration effect) caused by larger-sized SBCs ( Table 4).
In a previous clinical study, the authors advocated that the presence of SBCs leads to degeneration and wasting of cartilage and also increases the risk of surgical knee replacements (Tanamas et al., 2010). Our findings provide the biomechanical bases of how SBCs of different parameters (size, numbers, and locations) can cause degenerative effect on cartilage (i.e., increased stress effect). These findings also strengthen the hypothesis proposed by Radin and Rose (1986) that the progression of articular degeneration is associated with perturbations in underlying bone.
The maximum strain was noted on the cyst surface facing the joint line. Inner sections of the cysts showed that the maximum strain was at the top of the cyst (represented by " * " in Figure 5). As this part of the cyst was very near to articular cartilage, the most deteriorating effect of SBCs on cartilage was in this position. This intra-cystic focal effect in 12 mm-sized cysts (C-1 and C-3) was nearly 4 times greater than that of SBCs in sub-groups B-1 and B-3. Only a 1.4-fold increase in strain was produced by cysts at lower joint position in group B (B-2) than the corresponding cyst in group C (C-2). The strain spread was toward peripheral regions of the cysts (dotted circle in Figure 5), which means that cysts may expand toward their peripheral regions. SBCs located at the lower joint positions (sub-groups 2 of groups B and C) and SBCs located in posterior tibial halves demonstrated lower strain values than the same-sized SBCs at higher joint positions. This difference is clearly seen in Figures 8A,B, where the green line represents the strains in different groups. Sabokbar et al. (2000) implied that elevated pressure in bones surrounding the SBCs could evoke a macrophage response within subchondral marrow and contribute to bone resorption and cyst enlargement. Our results showed 1.6-fold increase in peri-cystic cancellous bone stress in group B than the non-cystic model (group A). This stress effect was increased up to 4.5-fold in SBC model of 12 mm size (group C). Increased intra-osseous stress may result in boney destruction and subsequent SBC expansion. The stress produced by cysts in cancellous bone in relation FIGURE 7 | Graphical representation relationship between cancellous bone stress (bone resorption) and strain.
Frontiers in Bioengineering and Biotechnology | www.frontiersin.org to strain is depicted in Figure 7. A similar phenomenon of cysts expansion has been previously described in the hip joint (Ondrouch, 1963;Durr et al., 2004). Ondrouch (1963) carried out photo-elasticity experiments and concluded that cysts in OA may develop as a result of stress-induced bone resorption. Durr et al. (2004) studied the etiology of SBCs in a 2D FE hip model and advocated that micro-fractures in subchondral bone may be the initial step in the formation of SBCs. Our results showed that the mean VMS in subchondral bone of SBCs was increased up to 4∼4.5-fold in sub-groups 1 and 3 of group B and nearly 16∼18.5fold higher in sub-groups 1 and 3 of group C as compared to the mean VMS in group A. In sub-group 2 of group B, only a Frontiers in Bioengineering and Biotechnology | www.frontiersin.org slight increase in mean stress was observed (p = 0.254). But a significant increase in mean stress was observed in sub-group 2 of group C (p < 0.001). Therefore, the location of SBC in relation to the joint line is also important in case of smallsized cyst to produce significant stress changes (p = 0.044 vs. p = 0.254). The blue lines in Figure 8B shows tremendous increase in subchondral bone stress in sub-groups C-1 and C-3 in contrast to the subchondral bone stress in sub-groups B-1 and B-3 in Figure 8A. This signifies the size effect of SBCs (12 vs. 4 mm). The different location of same-sized cysts (B-1 vs. B-2 and C-1 vs. C-2) also produced different stress effects as represented by red and blue lines for cancellous and subchondral bones, respectively (Figures 8A,B). We should also keep in mind that the repeated or increased loading beyond the body weight (singleleg standing) may further deteriorate the condition (further increase stress effect).

CONCLUSION
In conclusion, the bigger-sized SBCs in anterior and immediate position in the medial tibial cartilage produce increased focal stress effect in articular cartilage. Multiple cysts further deteriorate the condition by increased osseous stress effect and high tendency of peripheral cyst expansion.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the corresponding author on reasonable request.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Institutional Review Board of First Affiliated Hospital of China Medical University. The patients/participants provided their written informed consent to participate in this study.