WITHDRAWN: A Novel Robust Local Anisotropic Clustering Model for Tissue Segmentation and Bias Field Correction of Brain MR Image

<p class="0affiliationCxSpMiddle">This paper had been withdrawn by the authors, but has been inadvertently published by iJOE. It is hereby officially deleted from the issue.</p><p class="0affiliationCxSpLast">The International Journal of Online and Biomedical Engineering apologizes to the authors for the trouble caused.</p>


Introduction
Magnetic resonance (MR) image segmentation is an important step in MR image processing and analysis, and its results will directly affect diagnosis and treatment.Numerous models have been proposed, including clustering model, level set model, active contour model, and so on [1].The accurate segmentation results have important guiding significance for clinical diagnosis.However, misclassification is often caused by noise and bias field in brain MR image [2].Noise usually exists in MR image and it cannot be avoided during imaging process because of the impact of equipment.Bias field is generated by intensity inhomogeneity of the radio frequency field, which ap-pears as the smooth and slow change of the gray value of pixel for the same tissue.Hence, these defects in the process of brain MR image segmentation will directly lead to error segmentation [3].
Among all the segmentation methods, the most ripe and popular used is the clustering model that can be divided into hard clustering and soft clustering [4].In hard clustering, k-means clustering is a typical model of it.Soft clustering can be further classified into mixture models and fuzzy c-means (FCM), and the widely applied clustering method is FCM model.The advantage of FCM model is that it assumes pixels in the image pertain to more than one class, but it is noteworthy that FCM does not take any spatial information into consideration [5].
To enhance the robustness of FCM, many modified models incorporate spatial information into clustering.Ahmed et al. proposed the bias corrected FCM (BCFCM) model, which introduces the Euclidean distance between pixels in the neighborhood and its cluster center as the constraint term for the objective function of fuzzy clustering [6].However, BCFCM is not anisotropic, i.e., each pixel in the neighborhood has the same weight to the center pixel.Besides, the distance between neighborhood and cluster center must be recalculated during each iteration.
To decrease the time consuming and complexity, Klinids et al. proposed fuzzy local information c-means (FLICM) model, which introduces fuzzy parameters to directly manipulate the image and need not input a plethora of parameters [7].FLICM has strong robustness, but it always produces over-segmentation at the boundaries of different tissues and it does not take the bias field into consideration.
All above-mentioned are not consider that how to correct the bias field generated by intensity inhomogeneity.Hence, Li et al. designed the local intensity clustering (LIC) model [8] and multiplicative intrinsic component optimization (MICO) model [9], respectively.These two models applied the edge-based level set segmentation model and clustering model, respectively.Both LIC as well as MICO have strong ability to correct bias field existing in the image.MICO also takes full advantages of the property that bias field is smoothly and slowly changing.Although both models can effectively correct the bias field, they do not deal well with noise in MR image.Another popular model is bias correction embedded fuzzy c-means (BCEFCM) proposed by Feng et al. [10], which add the CUDA accelerated non-local means denoising model to reduce noisy pixels in MR image.They regard the original image as the product of real image and bias field and then convolve bias field with a kernel function in each iteration.
In this paper, we proposed a novel robust local anisotropic clustering (called RLAC) model to segment brain MR image and correct bias field simultaneously.First a weighting factor is constructed to make each pixel has anisotropic property, and each pixel in the original image is replaced by image patch, which can improve robustness to noise.Second, a multiplicative framework is used to correct the bias field, in which the observed image can be regard as a product of real image and bias field, and noise is viewed as an additive component.It is noteworthy that noise has been almost removed by the image patch in the weighting scheme.According to the property of the bias field, we represent it using a linear combination of the orthogonal Legendre polynomial, and it can ensure that the bias field can be correctly estimated.

2
Related Work

Fuzzy local information c-means
Krinidis et al. proposed the fuzzy local information c-means (FLICM).FLICM does not need to manually input a number of parameters and directly introduce fuzzy parameters to the image.The expression of the fuzzy parameter is given by (1) where represents the pixel within the neighborhood, represents the cluster center, represents the membership matrix, and represents the Euclidean distance between pixel and pixel.After introducing fuzzy parameters, the objective function of FLICM can be represented as (2) where the expressions of is defined by (3) the expressions of is given by ( FLICM has strong robustness to noise, but FLICM is not always effective in practical applications and incorrect segmentation may occur at the boundaries of different tissues.The reason for this situation is that the gray values of two adjacent tissues differ greatly, so a wrong value of will be obtained in both tissues.

Bias correction embedded fuzzy c-means
The FLICM model can only remove noise, but cannot correct the bias field.To accurately segment images and correct the bias field corrupted by both noise and intensity inhomogeneity, Feng et al. improved the FCM and then proposed the bias correction embedded fuzzy c-means (BCEFCM) model, in which the CUDA accelerated ki G ( ) non-local means denoising (NLMD) model is introduced to remove noisy pixels.The NLMD method was first achieved by Buades et al., which can preserve the original information effectively and will not blur the edges of the image.They do not use the neighbourhood information of the given pixel, but average the similar pixels in whole image domain, and the estimated values are given by (5) where is the noisy image.According to the similarity from to , They use that defined as the weights to replace original pixel with the weighted average of , which can be written in the following equation (6) where represents the filtering coefficient.
Owing to the computational time of NLMD, they incorporate it into the Compute Unified Device Architecture (CUDA).Finally, they introduce a series of positive weighting coefficients to segment the image into clusters.The energy minimization function is defined as (7) The BCEFCM model first calculates the weight coefficients using NLMD model, which is accelerated with CUDA.This model can accurately remove the noise from noisy MR image, but the bias field estimated by BCEFCM does not satisfy the smoothly and slowly varying property.

Robust Local Anisotropic Clustering
In some clustering models, they not pay any attention to the local information of image, and the results often sensitive to noise and outliers [11].Thus, RLAC model is proposed in this section to segment the brain MR image corrupted by noise and intensity inhomogeneity.As shown in Figure 1, RLAC model consists of two main parts: anisotropic weighting scheme and bias correction model.From the anisotropic weighting scheme in the first and second column of the Figure 1, it can be seen that the original image with severe noise and intensity inhomogeneity (Figure 1(a)) only in the foreground consisting of white matter (WM), grey matter (GM), and cerebrospinal fluid (CSF), thus the corresponding histograms of the original image not have the well-defined peaks.In contrast, noise removed image (Figure 1(b)) using the image patch is more noise free than original image, but it still has no well-separated peaks in the corresponding histograms because of the intensity inhomogeneity.Then, the bias x y y y x y y

Anisotropic weighting scheme
In the neighborhood of center pixel，each pixel has the different effect to center pixel, pixel is similar with the center pixel will has bigger weighting and vice versa.Thus, we use the weight of each pixel in the neighborhood and the corresponding pixel to update the current center pixel.
We consider the 8-neighborhoods for each pixel x.Assuming that the pixel within its neighborhood is r, then .We first calculate the gray mean square error between each pixel r and other pixels in the neighborhood.By means of the mean squared error, we can get the difference of gray value for all pixels in the neighborhood.The formula of mean square error is where and represent the gray value of r and the gray value of other pixels in the neighborhood, respectively.
shows the number of pixels within the neighborhood , i.e., .We then use the radial basis function to generate the weight for each pixel (9) Next, we normalize the weight of each pixel in the neighborhood (10) Through (10), we can divide the similar pixels in 8-neighborhoods into one class and update the gray value of center pixel according to this class that occupies the greater weighting.Finally, the gray value of center pixel is given by (11) The gray value is named the image patch, which used to replace each pixel in original image.

Bias field framework
In order to more accurately reflect the bias field and real image from the observed MR image, we assume that the observed image I is a multiplicative framework in a continuous domain.We can regard all the pixels in MR image as a product consisting of bias field, real MR image, and noise.The observed MR image is defined by (12) where is bias field.is the true image that reflects the intrinsic physical properties of different tissues, so under ideal conditions, pixel x in the uniform tissue should take the same value.
is the zero-mean additive Gaussian noise.According to the smooth and slow changing property of bias field, we express as the following form (13) where are the basis functions and are the optimal coefficients.In this paper, the linear combination of multiple order Legendre polynomial functions are used as the basis functions.Therefore, we can use a column vector to , M w w L describe the basis functions, which can be express as , where is the transpose operator.Then, the optimal coefficients also can be expressed as a column vector .Thus, can be rewritten as The

Energy formulation for RLAC
In this subsection, we propose an energy function based on the expression of noise removed image, bias field, and the true image to correct bias field and segment brain MR image, which can be defined by (16) where is denoising image after using image patch.We introduce the expressions of b and J described above for minimize this energy function.So the energy function can be rewritten as an equation about three independent variables, which are , , and , respectively.Thus, the energy equation can be expressed as (17) where is the weighting coefficient for cluster center.By minimizing the energy equation ( 17), we can obtain real brain MR image and bias field .Since exchange the order of integrals and summations does not change the equation and is a binary membership function in each .The energy equation can be rewritten as the following form (18) ( )

Energy minimization
Energy minimization is achieved through alternative iterations.The description of three variables in energy minimization of will be described separately.By minimizing the energy equation ( 18) with respect to and fixing and , The optimal solution of can be obtained by solving the partial derivative We can easily obtain the optimal solution of through the equation , which can be expressed as (20) Hence, the bias field can be obtained from the optimal solution of according to (14).
For fixed and , the optimal solution of can be obtained by the following equation (21) For fixed and , we use a binary function as the membership function of each tissue in this paper, which can ensure that the noisy free pixels can be correctly classified, so the function can be expressed as where .
In this paper, we set and four order orthogonal Legendre functions are used as the basis functions, which can be represented by , , , F w,c,u  w  c u 1, arg min ( ( )) ( ) 0, arg min ( ( )) , where and are directional components of image , and the number of terms of basis functions is 15.Then the specific implementation process of RLAC can be summarized as the pseudo code in Table 1.

Experimental Results
To show the robustness of RLAC and the ability to correct the bias field, we first use brain MR images only corrupted by severe noise to prove the robustness of RLAC.Then the RLAC model is used to correct the bias field on brain MR images only with intensity inhomogeneity.Finally, we compare RLAC with state-of-the-art models on brain MR images corrupted by both noise and bias field.All experimental results are qualitatively compared with state-of-the-art models, and also tested in Matlab R2016a on the computer with Core(TM) i5-7300 2.50GHz CPU, 4GB RAM, and Windows 7 operating system.

Robustness of RLAC
To demonstrate the robustness of RLAC to noise, the synthetic brain MR images with 9% noise and without inhomogeneity obtained from BrainWeb [12] are used as the experimental images and compared with BCFCM, FLICM, and MICO.It can be seen from the Figure 2 that the upper row shows the 67th slice image and the 96th slice image are showed on the lower row.The second column shows that noisy pixels are clearly present in the segmentation results of BCFCM since the local spatial information used by it is not anisotropic.FLICM utilizes the spatial Euclidean distance in the neighborhood, which will easily result in over-segmentation in adjacent different areas.MICO does not take the local spatial information into consider and hence sensitive to noise, so it will produce an incorrect segmentation for noisy brain MR tissue.Therefore, we can clearly observe that many outliers exist in MICO segmentation results.RLAC utilizes the weight of each pixel in the 8-neighborhoods to update the current center pixel, so it can identify the noisy pixels in the neighborhood, which can decrease the influence of noise.From the above experimental results, the robustness of RLAC is slightly superior to FLICM, and significantly superior to BCFCM and MICO.To more objectively and accurately compare the performance of the four models, the quantitative evaluation is used to analysis the segmentation results of the four models.We can also obtain ground truth (GT) from BrainWeb, which is the standard segmentation result of corresponding image.The Jaccard similarity (JS) coefficient can be used to compare the similarity of two sets.It is defined by (23) where and are two segmentation results that need to be compared.represents the number of pixels in corresponding area.We compare the segmentation results obtained by each method with GT.The higher JS is, the more accurate the seg- • mentation results.We calculate the JS coefficient for GM, WM, and CSF acquired from the BCFCM, FLICM, MICO, and RLAC model.Figure 3 shows the boxplots of the segmentation results.We can distinctly observe that RLAC has the highest JS compared with other three models and no outliers exist in the JS of RLAC, which shows the strong stability of RLAC.

Performance of RLAC in correcting bias field
In this subsection, the axial-sectioned, sagittal-sectioned, and coronal-sectioned brain MR images with severe intensity inhomogeneous are used to demonstrate the performance of RLAC in correcting bias field.From the Figure 4, the upper row shows the axial-sectioned image, the middle row shows the sagittal-sectioned image, and the lower row shows the coronal-sectioned image.It can be seen from the second column, the bias field satisfies the smoothly and slowly varying property since we use a multiplicative framework with a clearly expression of bias field, and also represent is using a linear combination of the four order Legendre polynomial.From bias field correction images in third column, we can observe that each tissue has the average intensity after using RLAC model.As shown in the last column, it can be observed that WM, GM, CSF, and background are distinct segmented into four classes.
We also compare the histograms of original image and bias field correction image to prove the improvement of image quality.As shown in the histograms of Figure 5, the histograms of original images for axial-sectioned, sagittal-sectioned, and coronalsectioned are shown in the first row, and the histograms of bias field correction images for axial-sectioned, sagittal-sectioned, and coronal-sectioned are shown in the second row.We can observe that original pixels intensity of GM, WM, and CSF are overlapped with each other which means that we are unlikely classify the original image into four classes accurately.There are no well-defined peaks that represent different tissues in original histograms because of severe intensity inhomogeneity, but there exist the obvious peaks in the histograms of bias field corrected images.From the above comparison, we can observe that RLAC can achieve the better performance on brain MR images with severe intensity inhomogeneous.

Effectiveness of RLAC
In this subsection, the axial-sectioned brain MR images with 7% noise and 80% bias field obtained from BrainWeb are used to verify the effectiveness of RLAC, and the comparison results with two state-of-the-art models are given in Figure 6.The segmentation results of RLAC, LIC, and MICO are displayed in the first, second, and third rows, respectively.The bias field obtained by LIC model does not satisfy the slowly and smoothly varying property.MICO model cannot effectively remove the noisy pixels existing in the axial-sectioned image.RLAC model uses a four order orthogonal basis functions to represent the bias field, so bias field satisfies the slowly and smoothly varying property and can also obtain an accurate segmentation results.In addition, we can intuitively observe that the noisy image after using RLAC can clearly distinguish the boundaries between different tissues.
The visual results do not fully demonstrate the performance of RLAC.We finally calculate the JS coefficient for GM, WM, and CSF obtained from LIC, MICO, and RLAC, respectively and given the results in Table 2.The noise ranges from 3% to 7% and the intensity inhomogeneity ranges from 40% to 80%.It is noteworthy that the better match between results and GT images will lead to higher values for JS.We can clearly observe that the proposed RLAC model gives rise to better JS coefficient than LIC and MICO, appearing robustness to different ranges of noise and intensity inho-

Conclusion
A novel robust local anisotropy clustering model proposed in this paper is used to segment tissue and correct bias field.We construct an anisotropic weighting scheme to improve robustness, which can make full use of local spatial information.In this scheme, original pixels are replaced by the image patch and also can preserve the details from noisy image.Then, we used a multiplicative framework to represent the observed image with true image and bias field, which can segment noisy brain MR image and correct the bias field simultaneously.A linear combination of the four order Legendre polynomial functions are used to represent bias field to ensure smoothly and slowly varying property.Energy function is then defined based on the expression of anisotropic weighting scheme, bias field, and real image.RLAC proposed in this paper can effectively remove the noisy pixels existing in the MR image and it has a higher performance than BCFCM, FLICM, and MICO.The quantitative evaluation compared with LIC and MICO has been proved the superiorities of RLAC for robustness and JS coefficients on numerous brain MR images.We also successfully applied RLAC to axial-sectioned, sagittal-sectioned, and coronal-sectioned images and obtained a high quality segmentation results.Furthermore, we can observe from the experimental results that RLAC can easily handle different degrees of bias field and noise.

Fig. 1 .
Fig. 1.The composition of RLAC model.(a) Original image and corresponding histogram, (b) noise removed image and corresponding histogram, (c) bias field corrected image and corresponding histogram, (d) bias field and segmentation result.

Fig. 5 .
Fig. 5.The histograms of original image and bias field correction image.(a) The histograms of axial-sectioned image, (b) The histograms of sagittal-sectioned image, (c) The histograms of coronal-sectioned image.

Table 1 .
The pseudo code of RLAC