K-means Clustering in Knee Cartilage Classification: Data from the OAI

Article historys: Received Sep 22, 2020 Revised May 29, 2020 Accepted Jun 29, 2020 Knee osteoarthritis is a degenerative joint disease which affects people mostly from elderly population. Knee cartilage segmentation is still a driving force in managing early symptoms of knee pain and its consequences of physical disability. However, manual delineation of the tissue of interest by single trained operator is very time consuming. This project utilized a fullyautomated segmentation that combined a series of image processing methods to process sagittal knee images. MRI scans undergo Bi-Bezier curve contrast enhancement which increase the distinctiveness of cartilage tissue. Bonecartilage complex is extracted with dilation of mask resulted from region growing at distal femoral bone. Later, the processed image is clustered with k = 2, into two groups, including coarse cartilage group and background. The thin layer of cartilage is successfully clustered with satisfactory accuracy of 0.987±0.004, sensitivity 0.685±0.065 of and specificity of 0.994±0.004. The results obtained are promising and potentially replace the manual labelling process of training set in convolutional neural network model.


INTRODUCTION
Osteoarthritis (OA) is the main culprit of chronic disabilities in the United States, 26.9 million of adults are estimated to be affected [1]. The statistic is believed to surge to 59 million by 2020 [2]. Meanwhile, 16% of the elderly population was estimated to be affected by knee OA [3]. Knee pain problem is becoming more common among adults aged 40 years and above. According to the Control of Rheumatic Diseases (COPCORD) study carried out in Malaysia, there were 23% of patients who aged over 55 years showed significant clinical symptoms while 39% among those who were over 65 years old [4].
Knee OA is a degenerative joint disease of whole knee joint in which all its articular cartilage structures are damaged. In early detection of the disease, the pathologic events are dynamic that matrix synthesis and repair will be increased while osteophytes start forming to stabilize the injured joint. Clinically, patients may start a series of rehabilitation programs or simply stop the activities which induce joint pain. Late in disease, OA is said to be a total joint failure that most of the joint structures have undergone irreversible pathologic mutation. The transition of knee OA from a dynamic to an irreversible pathologic changes differs greatly from people, in many persons, may never experience the disease too [5].
Human knee articular cartilage is a composition of dense extracellular matrix, which is made of water, type II collagen, proteoglycans, with numerous glycoproteins and other non-collagenous proteins [6]. Referring to Figure 1, both femur and tibia have thin cartilage layers that allow two bones to glide against each other essentially without friction. Progressive loss of articular cartilage causes the friction between the two bones to increase which then generates inflammation and triggers pain through the nerve endings in the joint space. Felson  used and bone density) and local biomechanical factors (muscle weakness, body mass index and joint laxity) contributing in increment of OA diseases [7]. Figure 1. Human knee anatomy. [24] Although there are several pharmaceutical treatments for OA patients who are facing severe pain and functional or mobility disability, it was still insufficient. Therefore, non-pharmacological interventions for instance, knee cartilage magnetic resonance imaging (MRI) segmentation is still a motivational force in managing early symptoms of pain and the consequences of physical disability [9] [20]. Generally, knee joint segmentation methods can be grouped into manual segmentation, semi-automatic segmentation or interactive segmentation, and finally fully automatic segmentation model. In this paper, we proposed a fully automatic knee cartilage segmentation framework with K-means clustering in which the coarse result is then adapted to Chan-Vese model for active contour to obtain the cartilage segmentation.

Related Work
Knee joint segmentation is adopted to conduct quantitative assessment of OA marker, for instance, bone deformation, cartilage thickness and volume and osteophytes formation. Among the existing radiography, MRI sequence is more suitable for quantitative assessment of cartilage status compared to other conventional radiography methods. MR imaging provides sufficient signal-to-noise ratio (SNR) and contrastto-noise ratio (CNR) with no significant artifacts such as geometric distortion and signal distortion [10]. Besides, high spatial resolution of MRI ensures sufficient pixels for cartilage thickness measurement such that articular cartilage is relatively thin, about 1.3 -2.5mm thick in a healthy knee [11].
Substantial effort has been contributed by the researchers around the world working on knee cartilage MRI segmentation for decades. The segmentation methods are different in their level of automation and the level of priori knowledge needed for the user to classify the region of interest (ROI) [8]. The methods which require high level of user interaction need low level of priori knowledge on the knee structure and vice versa. A good example of interactive segmentation model was presented to segment knee cartilage by using locally statistical level set method (LSLSM) and compute its thickness using normal distance [12]. Meanwhile, a low-level priori knowledge needed segmentation framework is proposed to segment the knee cartilage through canny edge detection to extract the edge of cartilage for further region of interest (ROI) masking [13]. The model requiring moderate priori knowledge and user interaction is best illustrated by Folkesson, who proposed to implement unsupervised k-nearest neighbours framework to cluster tissue by selecting specific features such as voxel position, raw and Gaussian smoothed intensities and intensity derivatives [15]. Yin described a new approach named Layered Optimal Graph Image Segmentation of Multiple Objects and Surfaces (LOGISMOS) that introduced multi-surface interaction constraints to inhibit oversegmentation of cartilages and bones [16]. A high priori knowledge of knee structure is needed in allocating the segmentation barrier for a good segmentation result. As can be seen in Figure 2, random walker interactive segmentation model reported in [6] required user interaction to place seeds to label background and manually label the tibial, femoral and patellar cartilages as foreground. The model aims to reduce the error caused by inter-and intra-observer variability and hence ensuring the segmentation result reproducibility.  Figure 2. User interaction on labelling a knee MR image.
A fully automated cartilage segmentation using multiple atlases built and local structural analysis method was proposed by June-Goo to determine the seed points location for graph-cut based method. The results are promising but there are drawbacks in computation time, 30 minutes for atlas building procedure and 10 minutes for fusion labelling and region adjustment procedures, and incomplete segmentation when the femoral cartilage and tibial cartilage are too closely juxtaposed [17]. Intervention of deep learning in evaluating knee MR images is reported in [14]. A depth of 5 convolutional encoder-decoder network is trained with substantial amount of labelled image by an experienced radiologist. The results obtained are convincing with trade-off of laborious label work and long model training time. Active shape model (ASM) overcomes the limitation of region-based segmentation methods such as active contour and region growing. ASM includes the knowledge of the shape of the ROI and forms a deformable shape model to fit into an interested region, for instance, extracting distal femur and proximal tibia bone for joint space width accessment in [21]. ASM initialization was obtained by placing landmarks according to the shape of the objects. The landmarks will then move to new points in the normal direction from their original points. The model will change its shape based on the eigenvalues and repeat the process until convergence [8]. Fripp et al. proposed 3D ASM in knee bone segmentation. The author highlighted the necessity of larger training databases so as to improve the segmentation results [22].

Problem Definition
Manual delineation of the tissue in interest by a single trained operator is very time consuming [18]. The segmentation performed in clinical routine, is irreproducible and its reconstruction times up to several hours [19]. In recent years, researchers are looking forward for semi-automated segmentation and fullyautomated segmentation methods to reduce the delineation time of knee cartilage. From the study, computeraided segmentation provides more consistent segmentation on ROI with less bias or variability compared to manual delineation of cartilage structure [19].
We utilize a fully-automated segmentation that combined a series of image processing methods. Initially, the knee MR images were pre-processed with contrast and intensity enhancement before smoothing filter and unsharp masking were applied. K-means clustering method was used to cluster the tissues, divide cartilage from noisy background including fat tissue, synovial fluid, bones and ligaments. In this study, sagittal knee MR images were tested with the proposed segmentation framework. Comparisons between manual segmentation and the proposed method with comparison of random walker segmentation model were discussed in section 3.

RESEARCH METHOD
All the MRI data are obtained from the Osteoarthritis Initiative (OAI). The MRI will be preprocessed to improve their contrast, reduce the existing noise and sharpen the edges.

MR Image Acquisition and Software Used
The knee joint MR images are provided by The Osteoarthritis Initiative. The images were captured using water excitation double echo steady-state (DESS) imaging protocol with sagittal slices at 3.0T. Given that imaging parameters for the sequence is TR/TE: 16.32/4.71ms, size: 384x384. MATLAB 2019a was used to do the image enhancement and segmentation.  Figure 3 shows the proposed framework used in this study. The DICOM format knee slice is converted into 16-bit grayscale TIFF file format that remains full image details. The image is smoothened with median filter of window size of 3 x 3 to remove most existing noises. Later, we apply unsharp masking to enhance the acutance of the boundaries of bones and cartilages.

Bone-cartilage Mask Generation
As a common problem in medical images, the appearance of MR image is dark in nature and the contrast of cartilage is indistinctive among the neighbouring tissues. To overcome the problem, Bi-Bezier curve contrast enhancement (BBCCE) [23] is conducted to enhance the global brightness while preserving the mean brightness of the image. Conventional histogram equalization is not prioritized in medical image processing as its mapping curve, also named as cumulative frequency density curve, experiences a sudden jump which could pull the intensity distribution naively that distorts the image quality. Therefore, with gentle nature offered by BBCCE, the knee MR image can be enhanced by preserving pertinent knee features while retaining its nature image appearance. Then, a seed is placed in femoral bone region. Obtaining the binary mask for the femoral bone, we dilate the mask with disk size of 30 to crop the bone-cartilage interface (BCI).

K-means classification framework
The bone-cartilage complex obtained from the previous stage is then input into the K-means classification framework. K-means clustering aims to partition the targeted observation into k cluster with the nearest mean. Given that a set of observations 1 , 2 , … , where every observation is a d-dimensional real vector. This clustering method would divide the n observations into k sets, which is less than n.
Each centroid defines its cluster. At first, each data point is assigned to its nearest centroid with squared Euclidean distance. For instance, if is the collection of centroids in set O, each observation x is assigned to their relative cluster based on: ( 1 ) Centroids, S will be recomputed by taking the mean of all observations assigned to the centroid's cluster: The clustering algorithm iterates until a stopping criterion is met, for instance, when there is no observation change clusters, the sum of the distances is minimized or the maximum number of iterations set initially is reached. In this study, k = 2 is defined to separate cartilage as the foreground while the other tissues and bones as background.
Morphological opening operation is carried out so as to disconnect the cartilage with the unwanted parts. Then, Chan-Vese for active contour is applied to recover the losing details that caused by the opening process. Finally, the cartilages are extracted which can be seen in Figure 3. The segmentation result will be tested with Dice Similarity Coefficient and its classification performance test.

Evaluation Metrics
In this work, 40 saggital knee MR images are tested with the proposed framework. Dice similarity coefficient (DSC) is used to verify the performance of our framework. As there is no ground truth provided by the database, manual segmentation of the images is normally done by an experienced operator. DSC refers to degree of agreement of both segmented results. Given that A is the segmentation result from the proposed framework while B is the manual segmented result by the operator: Sensitivity measures the performance of the proposed framework on classifying the cartilage pixels while specificity measures the ability of the framework to classify non-cartilage pixels. Accuracy measures the overall classification performance. Given that true positive (TP), false positive (FP), true negative (TN) and false negative (FN), sensitivity, specificity and accuracy can be defined as:

RESULTS AND ANALYSIS
Referring to Figure 4b, a BBCCE enhanced knee MR image can be seen. The cartilages are highlighted and appear to be more distinctive from the neighbouring fat tissues which could increase the accuracy of classification at the following stage. Clustering the bone-cartilage complex with k = 2, coarse cartilage information is extracted from the background. The segmented result can be contaminated with the fat tissues or synovial fluid which has homogeneous intensity with the hyaline layers. Therefore, the coarse result undergoes opening process to remove the unwanted parts and recover the lossed details with Chan-Vese active contour method with fixed of 200 iteration and the final segmentation can be seen in Figure 4c.
We compared our fully automated cartilage segmentation model with random walker interactive segmentation model [a] and the results from their performance evaluation are tabulated in Table 1 and Table  2. From 40 successful segmentations, the proposed method yields DSC of 0.689 ± 0.059, accuracy of 0.987 ± 0.004, sensitivity 0.685 ± 0.065 of and specificity of 0.994 ± 0.004. Meanwhile, random walker interactive model gives DSC of 0.694 ± 0.062, accuracy of 0.988 ± 0.003, sensitivity 0.678 ± 0.092 of and specificity of 0.994 ± 0.003. Notably that both frameworks acquire lower DSC value and sensitivity value indicating that both frameworks are sensitive to the intensity variation of the cartilages. As can be seen in Figure 5, hyaline cartilage will experience structure change with variant intensity distribution which can lead to wrong classification of cartilage pixels. In other words, the proposed model is sensitive and potential to detect cartilage lesion at early OA stage. Besides, both the segmentation models give high accuracy and specificity in cartilage classification.    As a trade-off of labelling work, the proposed framework classifies the bright cartilage from neighbouring dark tissues through intensity classification of k = 2. In Figure 6, misclassification of neighbouring tissue as patellar cartilage can be seen such that the intensity of the misclassified pixels is homogeneous with the cartilage pixels. To overcome this, morphological operations in breaking the connection and removal of unwanted region are suggested. Apart from that, a distance parameter can be added into the framework to compute the distance of the targeted pixels from the centroid. We can easily filter out the misclassified parts once it reaches the preset maximum distance.
In a nutshell, the proposed method that classifies the cartilages automatically is competitive with semi-automatic random walker segmentation model. The proposed model can be used to replace the tedious and time-consuming cartilage labelling work for deep learning model and active shape model which require large training database to ensure good classification results.

CONCLUSION
In short, 40 knee joint MR images from sagittal view are successfully segmented through the proposed framework to extract the cartilage information. We introduced a novel framework of combining gentle BBCCE contrast enhancement, bone-cartilage complex extraction with region growing and dilation of binary mask and finally classification of cartilage with K-means clustering and active contour implementation to recover the lost information. The proposed fully automated model is competitive with other semi-automated model and shows great potential in assisting the physicians in OA diagnosis and labelling cartilage for substantial amount of training set for deep learning classification model.