Association between subchondral bone structure and osteoarthritis histopathological grade

ABSTRACT Despite increasing evidence that subchondral bone contributes to osteoarthritis (OA) pathogenesis, little is known about local changes in bone structure compared to cartilage degeneration. This study linked structural adaptation of subchondral bone with histological OA grade. Twenty‐five osteochondral samples of macroscopically different degeneration were prepared from tibiae of 14 patients. Samples were scanned with micro‐computed tomography (μCT) and both conventional structural parameters and novel 3D parameters based on local patterns were analyzed from the subchondral plate and trabecular bone. Subsequently, samples were processed for histology and evaluated for OARSI grade. Each bone parameter and OARSI grade was compared to assess structural adaptation of bone with OA severity. In addition, thicknesses of cartilage, calcified cartilage, and subchondral plate were analyzed from histological sections and compared with subchondral bone plate thickness from μCT. With increasing OARSI grade, the subchondral plate became thicker along with decreased specific bone surface, while there was no change in tissue mineral density. Histological analysis showed that subchondral plate thickness from μCT also includes calcified cartilage. Entropy of local patterns increased with OA severity, reflecting higher tissue heterogeneity. In the trabecular compartment, bone volume fraction and both trabecular thickness and number increased with OARSI grade while trabecular separation and structure model index decreased. Also, elevation of local patterns became longitudinal in the subchondral plate and axial transverse in trabecular bone with increasing OARSI grade. This study demonstrates the possibility of radiological assessment of OA severity by structural analysis of bone. © 2016 The Authors. Journal of Orthopaedic Research Published by Wiley Periodicals, Inc. J Orthop Res 35:785–792, 2017.

In osteoarthritis (OA), subchondral bone becomes sclerotic and undergoes structural changes such as the formation of osteophytes and bone cysts. 1 It has been suggested that first tissue-level changes in OA would occur in the subchondral bone before any signs of degeneration in the overlying articular cartilage, 2,3 and thus contribute to the pathogenesis of OA. [4][5][6][7] In respective of increasing evidence of subchondral bone modifications with OA, these changes are generally not assessed until late OA stages with clinical diagnostics and histopathology.
Recently, it has been suggested that different OA phenotypes (post-traumatic, metabolic, ageing, genetic, pain) should be considered when studying OA and developing potential treatment or intervention. 8 Most of the animal models studied are post-traumatic and have shown a loss of bone in the trabecular compartment [9][10][11][12][13][14] and the subchondral plate, [12][13][14][15] whereas a thickening of the plate after trauma has also been reported. 9 The loss of the subchondral plate has been associated with increased thickness of calcified cartilage. 14 Most relevant models for ageing OA can be found, for instance, in Dunkin-Hartley guinea pigs and primates. These species typically have a thicker subchondral plate [16][17][18] with increased density and reduced porosity as well as increased trabecular bone volume fraction, anisotropy, plate-like shape, and transverse orientation. These skeletal changes are thought to precede cartilage degeneration. 18 In human (aging), OA a clinical hallmark is increased bone mineral density. 19,20 Histologically, it has been shown that in OA bone formation and volume are increased in subchondral trabecular bone with little resorption, 21 whereas presence of resorption pits and microtrabeculae in human OA samples have also been reported. 22 Changes in trabecular bone microstructure seem to follow a general pattern in various skeletal site: Tibia, [23][24][25][26][27] femoral head, 28 and neck. 29 These changes include increase in trabecular bone volume fraction, thickness and reduced separation, whereas the overall trabeculae shape becomes rather plate-like. [23][24][25][26][27][28][29] Moreover, altered trabecular orientation and increased anisotropy have often been reported, 29 though this effect has not always been obvious. 25,28 This controversy promotes the need for evaluating variations of local patterns in skeletal structures using alternative methods to better understand the bone changes in OA. 25,28,30 In previous studies, the subchondral plate in OA has been reported to be thicker and more porous. 22,31 To our knowledge, there is no study reporting the actual relationship between histological disease severity and subchondral bone plate modifications analyzed volumetrically using human sample material. Moreover, calcified cartilage has been shown to increase its thickness with OA process 32 but it becomes virtually absent in late OA. 22,32 Despite increasing evidence that the subchondral bone would contribute to the pathogenesis of OA, 4-7 most diagnostic methods are still focused mainly on cartilage degeneration and erosion, and skeletal changes are acknowledged only at late stages. Therefore, we hypothesized that subchondral bone changes occur in parallel with cartilage degeneration. To test this hypothesis, subchondral bone changes in both trabecular bone and subchondral plate were evaluated as a function of disease severity, evaluated by OARSI grading system. 33 Subchondral bone changes were analyzed by two methods: (1) traditional morphological analyses in whole tissue volume, and (2) analyses of local bone patterns using a volumetric local binary patterns (LBP)-based analysis, which is a promising method to study distribution, orientation, continuity, and complexity of local structures.

Human Osteochondral Samples
Tibial plateaus were collected from 14 late OA-patients (age 76 AE 9 years: 2 males and 12 females) treated with total knee arthroplasty at Oulu University Hospital, Finland. Sample collection and the study protocol was approved by the Ethical Committee of the Northern Ostrobothnia Hospital District, Oulu, Finland (191/2000). Altogether, 25 osteochondral samples were prepared from tibial plateaus from areas visually classified into three categories in terms of macroscopical degeneration of the articular cartilage: (1) intact cartilage (N ¼ 11), (2) moderate cartilage erosion (N ¼ 7), and (3) partly or fully exposed subchondral bone (N ¼ 7). Cylindrical samples with a diameter of 6 mm were prepared from all macroscopic visual grades.

Imaging and Evaluation of Structural Parameters
Samples were stored in phosphate-buffered saline (PBS) for micro-computed tomography (mCT) imaging. Osteochondral cylinders were scanned with a mCT device at isotropic 27.8 mm voxel size (Skyscan 1172, Bruker microCT, Kontich, Belgium). Images were reconstructed and analyzed with manufacturer provided software. The subchondral plate was manually segmented from the trabecular bone. Both trabecular and cortical plate compartments were analyzed separately. One sample did not have any trabecular bone and could not be included in the analyses. Images were prethresholded to 60 and remaining pixels were binarized with adaptive thresholding using the mean of minimum and maximum values within a 5 pixel radius window. After binarization, images were subjected to morphological analyses and the outcome variables are described in ref. 34 The Volumetric Local Binary Patterns Method In this study, a new alternative method, known as LBP 35 was used for analysis of local bone changes as a function of histological OA grade. Briefly, in the LBP method, the neighborhood of a center pixel is evaluated for occurrences of equal or higher gray level values than in the center pixel. For each center pixel, a specific local pattern is then determined based on the locations of these occurrences which are called markers. The LBP method was applied to every pixel with at least one neighbor with a higher value than a fixed threshold in order to capture only information about the bone and avoid non-relevant data from empty space.
LBP-based analysis was performed with a custom-made algorithm (MATLAB version R2013b, MathWorks, Natick, MA). The LBP method used here is an improved version of the one reported in the literature. 36 While the previous analysis was performed in 2D, the method used here was developed to perform 3D volumetric analyses. For this purpose, instead of fitting 8 neighbors on a circle within the studied slice, a total of 26 neighbors were fitted on a sphere ( Fig. S1) within multiple slices, the studied pixel being the center of the sphere. The gray-level value of the neighbors was calculated by interpolation of the voxels values surrounding them. Histograms of patterns distribution were obtained from the analysis and the entropy of the local patterns was evaluated as well as the amount of different local patterns obtained. The histograms were converted to obtain the occurrences of each specific pattern as a percentage of all the patterns recognized within a VOI. The entropy of local patterns was calculated as follows: where P i contains the count of a specific local pattern i occurring in the image. If an image contains only one local pattern, the entropy of the patterns within the image is zero. The entropy here describes the randomness of local patterns in the volume of interest. The elevation angle of each pattern was calculated using principal component analysis. The elevation is the angle along the distal-proximal axis, representing the alignment of patterns toward the cartilage surface. The patterns without consistent orientation were ignored in the analysis, the following cases were considered as such: Not enough markers (<3), too many markers (>24), or when markers were located a further distance from the "centroid" than non-markers. As the mCT imaging acquisition protocol did not allow maintaining the same reference orientation on the other plans, the azimuth angle was ignored in the analysis.
Both the entropy and the homogeneity of the elevation were calculated to describe spatial distribution of the elevation angles and to provide an indication on continuity of the structures. An angle-level co-occurrence matrix (ALCM) was first calculated with volumetric rotational invariance to assess the spatial distribution of the elevation angle. The homogeneity of the elevation was then calculated from the ALCM as such: where each element (i,j) in the ALCM is the sum of the occurrences that a voxel with elevation i has a voxel with elevation j in its vicinity (one voxel-size radius). Eventually, a high value of homogeneity of the elevation corresponds to data with consecutive voxels with similar elevation angles.

Histopathology
After mCT imaging, cylinders were formalin-fixed and decalcified in EDTA. Paraffin-embedded blocks were sectioned to 5 mm and stained with Safranin O and Masson's trichrome. OARSI grading 33

RESULTS
All knees had large internal variation in OARSI grades depending on the site. The OARSI grades varied from 0.9 to 6.3 (mean AE SD ¼ 3.5 AE 1.7). High precision was noted for both histological and mCTderived morphological parameters as indicated by CV % ranges 0.81-3.19% and 0.21-3.80%, respectively (Tables S1 and S2). Results from histological sections (Fig. S2)   mineralized structures from histological sections were compared against mCT-derived plate thickness, mCT gave slightly smaller values, and the correlation (R ¼ 0.98, p < 0.001) was stronger. The majority of the structural and texture bone parameters showed strong correlations with OARSI grades (Table 1 and Fig. 2). The subchondral plate thickness and fractal dimension correlated positively with OARSI grade (R ¼ 0.79 and R ¼ 0.83, respectively) with decreasing bone specific surface (R ¼ À0.71). The trabecular bone volume fraction increased with OARSI grade (R ¼ 0.78) with decreasing bone specific surface (R ¼ À0.79). Also the trabecular number (R ¼ 0.80) and thickness (R ¼ 0.78) increased, while their separation decreased (R ¼ À0.74) leading to increased connectivity of trabecular network as indicated by increasing connectivity density (R ¼ 0.73) and decreasing trabecular pattern factor (R ¼ À0.78) with disease severity. Furthermore, trabeculae shape became more plate-like as indicated by a negative correlation between structure model index and OARSI grade (R ¼ À0.70). Conventional morphological analysis showed only moderately lower anisotropy (R ¼ À0.47) with OARSI grade.
The novel 3D LPB analysis showed ( Table 2) that in both the subchondral plate and trabecular bone compartments, a higher OARSI grade corresponded to a higher number of different patterns (R ¼ 0.89 and  R ¼ 0.69, respectively) and an increase in the entropy of patterns (R ¼ 0.80 and R ¼ 0.83, respectively). Samples with high OARSI grades presented a stronger proportion of patterns oriented longitudinally in subchondral plate and transversally in trabecular bone (Fig. 3) than samples with lower OARSI grades. The homogeneity of angle elevation was also negatively correlated with OARSI grade in subchondral plate (R ¼ À0.80) as well as trabecular bone (R ¼ À0.79) while elevation entropy increased (R ¼ 0.79 and R ¼ 0.69, respectively) with OARSI grade. These results indicate discontinuity in the orientations within consecutive patterns, suggesting sharper angles in the architectures and eventually higher complexity of the structures in general.

DISCUSSION
It is known that subchondral bone becomes more sclerotic in OA, 2 but there is little evidence where the earliest microstructural changes in OA actually start -within cartilage or subchondral bone. Although this long standing question remains to be further elucidated, our results indicate that subchondral bone modifications are strongly coupled to articular cartilage degeneration already at the lowest OARSI grades. Morphometric parameters of subchondral plate and trabecular bone were evaluated separately due to the specific structure and functions leading to possibly different adaptation to OA in these compartments. 38 Our results indicate that OA-related alterations are present in both.
Our results showed higher trabecular bone volume fraction with OA due to higher trabecular number and thickness and lower trabecular separation. In this study, one of the calculated morphological parameters was SMI. This metric gives an estimation of the shape of the trabecular structures, from plate to rod-like. Here, we observed a negative correlation between SMI and OARSI grade, corresponding to a plate-like appearance of trabeculae in advanced OA. This finding was expected based on previous morphometric guinea pig model of spontaneous OA 17,18 and human studies. [23][24][25][26][27][28][29] As an adaptation to OA, trabecular bone may increase its stiffness and strength [39][40][41] by changing its structure to resist higher loading due to cartilage loss.
For both the subchondral plate and trabecular bone, the specific bone surface decreased with OA. This parameter corresponds to the ratio of bone surface by volume, suggesting that the bone surfaces become more uniform in trabecular bone and smoother in subchondral plate with increasing OARSI grade. In a recent study, 14 a large variation in subchondral plate thickness was observed within healthy and OA knees, with a tendency for the plate to become thicker with disease progression. As expected, the subchondral bone plate became thicker with higher OARSI grades and this relationship was surprisingly linear. To further validate    42 increasing the porosity of the subchondral plate 43 and its hydraulic conductance. 32 Also in post-traumatic OA models, accelerated bone remodeling has been observed as indicated by an increased number of osteoclasts and osteoid in subchondral bone, increasing its permeability (angiogenesis). 14 Nevertheless, our results showed the best correlation between thickness of mineralized structures and OARSI only when subchondral plate and calcified cartilage where accounted for. However, the role of tidemark in OA pathogenesis is still unclear and should be further considered in future studies. 44 These results highlight the potential of evaluating OA severity from mineralized structures. In longitudinal in vivo studies, this could be achieved either by following changes in calcified cartilage using ultrashort echo time magnetic resonance imaging 45 or evaluating subchondral plate and trabecular bone properties from high resolution clinical CT. The latter may become even more of an attractive option when combined with novel image processing tools such as 3D LBP.
In this study, the feasibility of applying a novel 3D LBP image analysis method related to local bone patterns in mCT scans was evaluated. By using 3D LBP-based method, it was possible to derive information on distribution, orientation, and complexity of patterns. The primary advantage of LBP is its low sensitivity for monotonic grayscale variations. 46 This is important because LBP may be ideally suited for clinical CT or extremity CT images where beam hardening and partial volume effect can reduce the variation in gray-levels to an extent that individual structures cannot be accurately segmented using global or even adaptive thresholds.
To date, the LBP method has been applied solely in 2D to analyze local bone structures in OA. 36,47 As noted earlier, entropy of local patterns is a measure of the randomness of the different bone patterns appearing in the region or volume of interest. Recently, it has been shown that the entropy of local patterns was related to OA severity assessed via K-L grading scale, 47 which parallels our findings. Also LBP-derived structural local bone properties identified in this study are in accordance with previous fractal-based studies showing significant differences between controls and patients with OA. 48 The current study showed that the elevation of patterns with increasing disease severity indicated that the subchondral bone plate and trabecular bone have different characteristics, which is expected since these compartments have different morphology, physiology and mechanical properties. 4 The higher complexity and discontinuity of the local consecutive patterns (in terms of orientations) in both compartments is directly related to altered bone structure as more connections between the bony structures are created with advancing OA. 29 In the subchondral bone plate, the appearance of non-uniform patterns could indicate altered porosity and tissue heterogeneity. In trabecular bone, the degree of anisotropy was lower with higher OARSI grades as previously suggested. 25,28,29 The novel 3D LBP method clearly demonstrated higher entropy and a rather transverse orientation (Fig. 3). This altered orientation fits well with the early proposal by Radin and Rose about adaptation of subchondral bone to shear stress. 49 This study highlights the importance of subchondral bone changes in human OA. Most importantly, these results indicate that changes in articular cartilage and subchondral bone are likely to be coupled, at least in this sample set. Therefore, based on this data, the cartilage-driven initiation of OA cannot be supported. However, it should be mentioned that according to current understanding, OA is not just a single disease but instead a group of diseases with different phenotypes. 50 It may be that the current samples represents an ageing phenotype, where both cartilage and bone changes occur with relatively the same speed. On the other hand, post-traumatic OA, for example, could be a phenotype starting from the articular cartilage defect, leading to subsequent subchondral bone changes. However, at this stage, this can be only speculated and more research is needed to understand the etiopathology of OA and its different phenotypes.
A major strength of this exploratory study pertains in comprehensive analysis of bone and cartilage using histological methods as well as mCT analyses with both conventional and novel LBP-based methods. Furthermore, as this study focused on smaller volumes directly under histologically graded area, the link between bone morphometric parameters and OARSI grade was highly linear. A major limitation, however, remains concerning the tissue samples from late-OA patients scheduled for knee surgery. Given their nature, it cannot be known for sure whether the samples with low OARSI grades are truly representative of early OA bone and cartilage. To fully confirm this, a larger size of human samples from both symptomatic and asymptomatic donors is required. In summary, results obtained from volumetric standard bone parameters, as well as parameters from local analysis of bone patterns, were complementary. The high correlation between the OARSI grade and the subchondral plate thickness and the trabecular bone volume fraction indicate that bone sclerosis and cartilage degeneration are coupled. Furthermore, the novel LBP-based method is a potential tool for evaluating both quantitative and qualitative changes from low-resolution clinical images. Altogether, this study demonstrated the possibility to diagnose severity of OA using bone structural analysis, which eventually could be assessed in clinical settings.

AUTHORS' CONTRIBUTIONS
MF, JT, MTN, PL, SS contributed to study design. MAJF, OA, VT, JR, SK, MV, PL contributed to data collection. MAJF, JT, SK, KP contributed to data analysis. MAJF, JT, KP, PL, SS contributed to data interpretation. MAJF, JT contributed to manuscript preparation. All authors contributed to revising and approval of manuscript content. MF takes responsibility for the integrity of the work.