Efficient Johnson-SB Mixture Model for Segmentation of CT Liver Image

To overcome the problem that the traditional Gaussian mixture model (GMM) cannot well describe the skewness distribution of the gray-level histogram of a liver CT slice, we propose a novel segmentation method for liver CT images by introducing the Johnson-SB mixture model (JSBMM). The Johnson-SB model not only has a flexible asymmetrical distribution but also covers a variety of other distributions as well. In this article, the parameter optimization formulas for JSBMM were derived by employing the expectation-maximization (EM) algorithm and maximum likelihood. The implementation process of the JSBMM-based segmentation algorithm is provided in detail. To make better use of the skewness of Johnson-SB and improve the segmentation accuracy, we devise an idea to divide the histogram into two parts and calculate the segmentation threshold for each part, respectively, which is called JSBMM-TDH. By analyzing and comparing the segmentation thresholds with different cluster numbers, it is illustrated that the segmentation threshold of JSBMM-TDH will tend to be stable with the increasing of cluster number, while that of GMM is sensitive to different cluster numbers. The proposed JSBMM-TDH is applied to segment four randomly obtained abdominal CT image sequences, and the segmentation results and robustness have been compared between JSBMM-TDH and GMM. It is verified that JSBMM-TDH has preferable segmentation results and better robustness than GMM for the segmentation of liver CT images.


Introduction
Liver cancer is a common malignant neoplasm worldwide, and the incidence of primary liver cancer is still on the rise at the global level [1]. e accurate understanding of the shape of the liver, the location and size of lesions in the liver tissue, and the relationship between the liver and surrounding blood vessels can help doctors to develop more effective treatment options. In addition, accurate liver segmentation is also conducive to three-dimensional reconstruction for liver and virtual surgery. At present, manual delineation of each slice by experts is still the standard clinical practice for liver demarcation [2]. Because the segmentation of organs and lesions has to be carried out layer by layer in CT slices, it is pretty cumbersome and time-consuming for doctors or experts to do this repetitive work.
In practice, CT slices are grayscale images. e pixels in the image reflect the X-ray absorption coefficient of the corresponding voxels. Black areas represent low absorption areas, that is, low-density tissues and organs of human body, such as lungs; white areas represent high-absorption areas, that is, high-density body parts, such as bones.
e absorption value in CT images is given in Hounsfield units (HUs). Compared with ordinary X-ray images, CT images have a higher density resolution. erefore, CT images can better show organs composed of soft tissue, such as lung, liver, gallbladder, pancreas, and pelvis, and can well demonstrate pathological changes in the tomography image. At present, the original pixel size of the abdominal CT image is 512 × 512, and the HU value of the liver varies widely. For example, a healthy liver has smooth contours and uniform density, and the absorption values are 60 ± 6 HU (or 64 ± 5 HU). However, when the liver has about 80% steatosis, the absorption values will be reduced to about − 50 HU. In contrast, the density of the liver parenchyma will increase due to the accumulation of iron in patients with hemochromatosis, whose CT scan showed that the liver parenchyma was clearly dense and bright, with an absorption value as high as +140 HU (so-called white liver) [3]. ere is no doubt that it is quite difficult to accurately segment the liver within such a large gray range.
To solve the above problems, many scholars have proposed a variety of liver CT image segmentation methods [4][5][6][7]. Some of these algorithms require human-computer interaction, such as active contour [8] and Livewire (intelligent scissors) [9,10]. Some are semi-automatic, for example, graph-cut [11] and region growth [12,13] methods. Some other methods focusing on fully automatic segmentation include statistical shape models (SSMs) [14,15] and thresholding algorithm [16]. Neural network algorithms aiming to achieve automatic feature extraction have also been applied to the segmentation of medical images in recent years [17,18]. e ultimate common goal of different methods is to segment images accurately and automatically, but this goal is still a bottleneck problem in liver CT image segmentation, due to the complexity of abdominal CT images and the differences between different liver morphologies. In this study, we aim to study on the finite mixture model (FMM), one kind of threshold segmentation algorithm, to improve the segmentation accuracy of liver CT images, and try hard to segment automatically at the same time.
According to the idea of the threshold segmentation algorithm, if the grayscale threshold of the liver in the CT slice can be accurately determined, it is possible to realize the liver segmentation automatically. In 1893, Karl Pearson made an experiment using the method of moments to fit a mixture of two normal components to the crabs' data, which proved the FMM could improve the accuracy of clustering [19]. Since then, FMM was adopted to improve the accuracy of threshold segmentation methods. In 1972, Chow and Kaneko applied FMM in medical images to segment the left ventricle from cine angiograms with two Gaussian distributions [20].
ere are two core points in FMM, one is the selection of the probability density function of the mixed components and the other is the parameter estimation of the mixture model. e most common mixed component probability distribution used in FMM is the Gaussian distribution, because in many cases there is a normal distribution in univariate and multivariate data. erefore, the Gaussian mixture model (GMM) has been widely used in the segmentation of the images [21][22][23][24]. In addition to the Gaussian distribution, gamma distribution, Student's t distribution, exponential distribution, and Rayleigh distribution also commonly appear in FMM [25][26][27], and their probability density diagrams with different parameters are shown in ). Although the peak shape of gray-level histogram has symmetry to some extent, its asymmetry is also very obvious, which is different from the symmetric characteristic of Gaussian distribution and the asymmetric characteristic of the exponential, Rayleigh distribution, etc. erefore, it is very difficult to fit the peak shape accurately with any single distribution shown in Figure 1.
To solve the problem, many researchers have focused on forming a mixture model using distribution functions that can better fit the shape of a single peak. e research ideas mainly concentrate on the following three kinds of mixture models.
e first kind of mixture model is forming an asymmetric generalized Gaussian distribution (AGGD) by introducing shape parameters or functions into the generalized normal distribution to describe the skew characteristic [28,29], so that it can describe not only a symmetric distribution but also an asymmetric distribution. However, the expression of the AGGD is complicated by embedding the gamma function. In the case of using the EM algorithm for maximum-likelihood estimation (MLE), all AGGD parameters in the mixture model are represented by highly nonlinear equations, which makes the numerical solutions cumbersome and sensitive to initial EM values [30].
Combining Gauss with other distributions to form a new model is the second idea. Wilson selected two Gaussian distributions and one uniform distribution to fit the lowgray and high-gray regions of the brain MRA data histogram, respectively [31]. Hassouna proposed a linear combination of a finite mixture model using one Rayleigh distribution and two Gaussian distributions [32]. Hence, for different problems, people need to determine in advance which existing models can be used to form a new probability distribution model, and the number of each distribution also needs to be determined in advance. erefore, this kind of model is not flexible, and it is also difficult to realize segmentation automatically. e third way is to combine the components of the mixed model with the non-Gaussian distribution. Lee and McLachlan introduced a finite mixture of canonical fundamental deviation t (CFUST) distributions for asymmetric and possibly long-tailed clusters [33]. Sefidpour and Bouguila proposed and investigated the segmentation of spatial color images using the Dirichlet and Beta-Liouville distributions [34]. e normal inverse Gaussian distribution (NIG) is chosen by Karlis and Santourian to deal with skewed subpopulations [35]. Franczak et al. studied the asymmetric Laplace distribution (ALD) for clustering and classification [36]. NIG and ALD belong to the family of generalized hyperbolic (GH) distributions designed by Barndorff-Nielsen [37]. Browne and McNicholas extended a special case for the generalized hyperbolic distribution [38]. Wraith and Forbes studied the properties of these distributions in multiscale and their application in multivariate clustering [39]. Although there are many combination methods of non-Gaussian distribution, only a few methods are used for the segmentation in medical images, especially in liver CT images. erefore, the study on the segmentation of liver CT images still needs further research.
Considering the asymmetric skew characteristic of the gray histogram of the liver CT image, and inspired by the skew characteristic of Johnson-S B distribution, we propose a novel mixture model with the combination of Johnson-S B distribution to segment the liver CT image. In the following paper, the Johnson-S B mixture model and the optimized parameters with the EM algorithm are introduced firstly in Section 2. Secondly, Section 3 gives the implementation details of the segmentation algorithm based on the Johnson-S B mixture model, and the effects of cluster number on the segmentation threshold are also analyzed and compared between GMM and Johnson-S B mixture model in this section. en, the segmentation experimental results of four randomly obtained abdominal CT image sequences from different image databases are provided in Section 4. Finally, conclusions are drawn in Section 5.

Finite Johnson-S B Mixture Model
In 1949, Johnson deduced a curve system called the Johnson system, which contains Johnson-S B , Johnson-S L , and Johnson-S U . e symbol S L means "log-normal system," S B means "bounded system," and S U means "unbounded system" [40].
e Johnson system can closely approximate many continuous distributions by one of the three functional forms, so it is very flexible to fit variety curves. Many of the commonly used continuous distributions, such as normal, log-normal, gamma, beta, and exponential, are special cases of the Johnson system; therefore, it has more advantages to fit curves with the Johnson system than any other single distribution [41].

Johnson-S B Distribution.
Johnson-S B is one of the three distributions of the Johnson system. It corresponds to the distribution of a continuous random variable x in which a particular transformation is applied to obtain a normal distribution. e transformation is as follows: where x is a given continuous random variable. In this study, x is the grayscale value of the pixel in the CT image. x ∈ (ξ, ξ + λ), ξ � min (x), λ � max (x) − min (x), and c and δ are shape parameters, δ > 0, c ∈ (− ∞, +∞). Z is a standard normal random variable, and its probability density function is as follows: We write where z is the inverse function of y.
According to the transformations of continuous random variable, en, the probability density function (PDF) with regard to x is as follows: Figure 3 shows some typical probability density function curves of Johnson-S B with different parameters. e range of horizontal coordinates is from 0 to 255, and the vertical coordinates are the corresponding probability density function values. Figure 3(a) presents the curves with different δ and c � 0, and it is a normal distribution in the middle. Figures 3(b)∼3(e) present the curves with different δ and c � 0, respectively, and they have better skew characteristics, especially at the two ends of the abscissa. It can be seen that c controls the position of the function, and the distributions of the function are normal distribution, negatively skewed distribution, and positively skewed distribution when c � 0, c > 0, and c < 0, respectively. e finite mixture model (FMM) refers to the linear superposition of distribution functions of the same type but with different parameters. In the discrete case, the probability density function of a finite mixture distribution can be expressed as a p-dimensional random vector X [42].
where x ⊂ X is the pixel grayscale value in the CT image, (3) defines a density that is called a K-component finite mixture density. Here, the weight of the kth mixing proportion Φ k satisfies the following relations: and According to Eqn (3), the probability density of a Johnson-S B mixture model is defined as follows: Here, is the probability density of a random variable X for a Johnson-S B distribution with the parameters δ k and c k .

Optimizing Parameters with EM
Algorithm. e vector of parameters Θ typically introduced by the loglikelihood function is defined as follows: Here, x i is the ith discrete grayscale value and N is the number of discrete dots of CT image. e detailed derivation of the parameters given in equations (6)-(8) is presented in Appendix A.

Implementation Details.
ere are mainly four steps to implement the segmentation algorithm based on J SB MM, and the details of each step in Algorithm 2 are introduced as follows.

Obtaining the Approximate Gray Value (LV_A).
e Hounsfield unit (HU) value of liver tissue varies from patient to patient. In addition, the X-ray tube of the CT machine will age with longtime use, which results in the intensity decrease in the X-ray source. erefore, the liver HU values of different CT image sequences are usually different. Clearly, it is not possible to identify all liver scan sequences with the same HU range value.
We can make full use of the continuity of CTscan section images to solve this problem. CT scans require patients to hold their breath during data collection, and this process usually takes no more than 20 seconds [43]. erefore, in the same CT image sequence, the range of gray values between all liver slices varied a little. It is feasible to obtain the liver HU value of this sequence using the slice image with the largest liver area.
calculate the possibility of each xi coming from the kth submodel, based on the current parameters θ k .
(4) M-step: optimize the model parameters Φ k , δ k , c k of the new iteration by using maximum likelihood.
(1) Obtain an approximate grayscale value of the liver (LV_A) of a given CT image sequence. images of the liver CT scan sequence, respectively. It can be seen that the liver is almost completely surrounded by the ribs, and the liver image at about 1/3 of the entire sequence is the largest one [44]. For example, if a complete sequence of CT scans of the liver has 90 slices, the liver slice near the 30th slice will be more significant than others.

Finding the Maximum Gray Value (LV_M).
Although liver HU values do not change significantly during the scan process, there also exist some slight changes. To segment each slice accurately, it is necessary to further locate the maximum grayscale value of the liver (LV_M) close to the LV_A on the gray-level histogram of each slice.
For instance, Figure 6(a) shows the slice of the upper part of the liver CT sequence shown in Figure 5(c). Figure 6(c) displays the slice of the lower part in the same sequence. Figure 6(b) gives the gray-level histogram of Figure 6(a), and the gray-level histogram of Figure 6(c) is shown in Figure 6(d). Because the value of LV_A shown in Figure 5 Here, it needs to be stated that we should choose the peak value closest to 178 as LV_M, because 178 is determined in Figure 5(d).

Determining the Segmentation Points.
e position of peak and skewness characteristics of the gray-level histogram of liver slices are uncertain and nonuniqueness, so it is difficult to take full advantage of the J SB MM in describing skewness if J SB MM is applied directly to fit the whole graylevel histogram. To make full use of the advantages of Johnson-S B in describing the skew characteristics, we divided the gray-level histogram into the left and right parts at LV_M to ensure that the gray level of the liver is just at the boundary of the histogram, to improve the Johnson-S B segmentation accuracy. e LV_M of Figure 6(c) is 176. At the LV_M point, the gray-level histogram shown in Figure 6(d) is divided into the left part and the right part, and the grayscale values of the segmentation points are obtained using J SB MM to fit the gray-level histograms. Figures 7(a) and 7(b) give the fitting results of the left part and the right part, respectively. e grayscale value of the segmentation point of the left part is 157, which is the intersection of the two rightmost fitting curves shown in Figure 7(a). e grayscale value of the segmentation point of the right part is 190, which is the intersection of the two leftmost fitting curves shown in Figure 7(b). However, if there is a peak of one fitting curve appearing between the intersection point of the other two fitting curves and the LV_M, the grayscale value of this peak will be taken as the segmentation point, which can get better segmentation results. For example, the grayscale value of the segmentation point of the right part should be 188, because the grayscale value of the peak is 188, which appears before the intersection of 190. erefore, the range of the grayscale value for liver segmentation in Figure 6(c) is set from 157 to 188, which is the range of the segmentation thresholds for Figure 6(c).

Binarizing the Images and Processing with Mathematical Morphology.
Mathematical morphology is a tool for extracting image components. Erosion, dilation, and filling are three basic morphological set transformations. ese transformations involve the interaction between image and structuring element [45]. Figure 7(c) is obtained by binarizing Figure 6(c). e final boundary of the liver segmentation can be drawn as shown in Figure 7(d) after the procession of filling, erosion, and expansion algorithms in mathematical morphology.

Effect of Cluster Number on the Segmentation resholds.
A key point of FMM is how to select the cluster numbers to realize a better fitting, and underfitting or overfitting may occur if the cluster number is selected inappropriately [46]. In this study, the cluster number means the number of curves in a cluster used to fit the gray-level histogram. To analyze the effect of different cluster numbers on segmentation thresholds, we choose two different slices to be segmented by GMM and J SB MM in the following Section 3.2.1 and Section 3.2.2, respectively. e initial location of the grayscale value is another consideration to affect the segmentation threshold. In this study, the initial locations are evenly arranged within the range of the grayscale values of the entire image.

Segmentation
resholds of GMM with Different Cluster Numbers. Table 1 gives the segmentation thresholds of the upper, middle, and lower slices from two different CT image sequences. By analyzing the segmentation thresholds with different cluster numbers n, it can be found that the segmentation thresholds are changeable with the change in n, while they do not tend to be stable with the increasing n. In this study, n denotes the cluster numbers used to fit the gray-level histogram, and n is the same as K in Algorithm 1.
To illustrate the effect on the segmentation results with different cluster numbers, we give some segmentation results of the upper slices of the sequences S1 and S2. Here, Figure 8 shows the threshold segmentation results of the upper slices with n � 6, 8, 10, 12, and 14 using GMM, respectively. e first row and the second row give the segmentation results and thresholds for the upper slice of the first sequence S1. e third row and the fourth row give the results and thresholds for the upper slice of the second sequence S2.
By comparing the segmentation results in the first row, it is not difficult to find that the results of the first three columns with n � 6, 8, and 10 are similar, because the three segmentation results include some other organs that do not belong to the liver. Although the segmentation result in the fourth column with n � 12 only includes the liver, the segmentation result is a little smaller liver than that of the fifth column with n � 14. erefore, the segmentation result with n � 14 is relatively more reasonable by comparing with the other four results for the first sequence.
Comparing the results of row 3 for the second sequence, it is obvious that the results of the first two columns with n � 6 and 8 are not correct, and the segmentation result of the third column with n � 10 is a little smoother than that of the fourth one with n � 12, while the segmentation result of the fifth column with n � 14 turns to be rougher boundary due to overfitting. us, the better segmentation result can be obtained when n is 10 for the second sequence. As can be seen from Figure 8, the segmentation results of GMM are sensitive to cluster numbers, and it is not easy for GMM to use fixed cluster number to get accurate segmentation results, which is a difficulty in the application of GMM, especially for a large number of images that need to be segmented.
Here, it is needed to be stated that the red contours in the figures in Figure 8 indicate the error zones of the segmentation. It is the same meaning in the following figures. Figures 9∼11 give the segmentation thresholds and the results of the upper, middle, and lower slices of the first liver CT sequence S1, and here, the segmentation thresholds are obtained with J SB MM-TDH, which means the gray-level histogram is divided into left and right parts, and it is denoted as thresholds of dividing histogram (TDH). While thresholds of whole histogram (TWH) mean that the segmentation thresholds are obtained with J SB MM under the whole gray-level histogram, and it is denoted as J SB MM-TWH. e segmentation results of J SB MM-TWH are not given in Figures 9∼11, and only the segmentation thresholds are summarized in the TWH column given in Table 2 e segmentation thresholds with n � 6∼15 for the upper, middle, and lower slices are summarized in the corresponding TDH column given in Table 2. It is obvious that the segmentation thresholds almost tend to be stable when n is bigger than 12. While the segmentation thresholds in TWH column are also changeable, even n is bigger than 12, which is similar to the trend of GMM. By analyzing and comparing the range of segmentation thresholds between J SB MM-TDH and J SB MM-TWH with the same n, it can be seen that most of the upper limit values are close to each other, but the lower limit values of J SB MM-TDH are smaller than that of J SB MM-TWH, and the differences are mostly around 10.   Analyzing and comparing the segmentation results of (a)∼(d) shown in Figures 9∼11, respectively, the better segmentation results can be obtained for the upper, middle, and lower slices when n is 13. erefore, considering the trend of segmentation thresholds given in Table 2, we determine to take n � 13 as the fixed cluster number for J SB MM-TDH to segment different CT image sequences in the following segmentation experiments in Section 4.

Segmentation resholds of J SB MM with Different Cluster Numbers.
As a matter of convenience in comparison, Table 2 also gives the segmentation thresholds of GMM for liver CT sequence S1 shown in Table 1. By analyzing and comparing the range of segmentation thresholds between J SB MM-TWH and GMM with the same n, it is found that most of the upper and lower limit values are close to each other, which means the segmentation results of J SB MM-TWH are similar to that of GMM. erefore, the segmentation results of J SB MM-TDH are both better than the results of J SB MM-TWH and GMM. Table 3 gives the segmentation calculation time with different cluster numbers as shown in Figures 9-11. e segmentation time is less affected by the different cluster numbers, and it almost tends to decrease for the middle slice of Figure 10 with a large liver area. In this study, only the parameters Φ, δ, and c are used to calculate the threshold value of analytic solution, so the segmentation speed is quite fast, and the effect of segmentation time with different cluster numbers can be ignored.

Segmentation Experimental Results
e J SB MM-TDH was applied to segment four randomly obtained abdominal CT image sequences from different image databases [47,48], and the segmentation thresholds and results of J SB MM-TDH are compared with that of GMM, as shown in Figures 12∼14, which are the upper, middle, and lower slices of the four different CT image sequences, respectively. In these figures, the first row shows the original image from four different CT image sequences, the second row gives the threshold segmentation results of GMM, and the third row is the binary images obtained according to the segmentation results of GMM. e fourth row shows the binary images obtained depending on the segmentation results of J SB MM-TDH. e fifth row is the segmentation results of GMM, and the sixth row gives the results of J SB MM-TDH. e seventh row and the last row provide the left and right segmentation thresholds of J SB MM-TDH, respectively. Table 4 gives the comparison of segmentation thresholds and quantitative evaluation for Figures 12-14, respectively. Here, it needs to be stated that the segmentation results of GMM are fitted with different cluster numbers, which is the best one chosen from the segmentation results with n � 6∼15 for each slice for the sequences S2, S3, S4, and S5, respectively, while all segmentation results of J SB MM-TDH are fitted with the same cluster number n � 13. However, most segmentation results of J SB MM-TDH are better than those of GMM, and the other few results are similar to each other. On the whole, it is illustrated that not only the segmentation results of J SB MM-TDH are better than those of GMM, but also the cluster numbers can be fixed at n � 13 without worrying about overfitting.
is is good for the implementation of automatic segmentation of liver CT images.
e Jaccard index and Dice coefficient are common indexes for quantitative evaluation of image segmentation.
e Jaccard index is a statistic used for comparing the similarity of sample sets. e Dice coefficient is another similarity measure index. e Jaccard index and Dice coefficient are calculated to quantitatively evaluate the segmentation results of GMM and J SB MM-TDH, respectively. e values of the Jaccard index and Dice coefficient are given in Table 4, and the values in the column of difference mean the difference value between J SB MM-TDH and GMM, that is, the Jaccard index (or Dice coefficient) of J SB MM-TDH minus that of GMM. e last row gives the average value. e maximum value of difference column in the Jaccard index is 0.1987, and the average value of this column is 0.0691. e maximum value of difference column in the Dice coefficient is 0.1863, and the average value of this column is 0.048. By comprehensive quantitative comparison, it is found that the segmentation results of J SB MM-TDH are better than that of GMM.
It is obvious that the segmentation thresholds almost tend to be stable when n is bigger than 12, which is shown as the bold values shown in TDH columns.     e website addresses of the data sets S1, S2, S3, S4, and S5 are provided in Appendix B.

Conclusions
e J SB MM-TDH with flexibly skewed characteristics proposed in this study is suitable for fitting the skewness distribution of the gray-level histogram of liver CT images. e parameter optimization algorithm employing EM and the implementation process of the segmentation algorithm has been given in detail. e effects of cluster number on segmentation threshold were discussed and compared for GMM, J SB MM-TWH, and J SB MM-TDH, respectively. It is shown that the J SB MM-TDH threshold will tend to be stable at cluster number 13, while the threshold of GMM and the threshold of J SB MM-TWH are similar and sensitive to different cluster numbers. e proposed J SB MM-TDH with cluster number 13 is applied to segment four random CT image sequences, and the segmentation results are compared with those of GMM. Analyzing the segmentation results and quantitative evaluations, it is further illustrated that J SB MM-TDH does not have the overfitting phenomenon with the increase in cluster number, which verifies that J SB MM-TDH has preferable segmentation results and better robustness than GMM. J SB MM-TDH makes it possible to realize the automatic segmentation of live CT image due to the robustness with the fixed cluster number. e J SB MM-TDH can be used not only for liver CT image segmentation, but also for other CT image segmentation as well.

Appendix
A M-step updates equations (6)∼(8) for solving the JSBMM parameters presented in Section 2.2.2, which are derived by maximizing the complete data log-likelihood Q with respect to each model parameter as follows: en, Let ln en, Φ k is solved by applying the method of Lagrange multipliers as follows: (A.18) Journal of Healthcare Engineering 19