simi-lar principle, the priori information from probabilistic atlases was used to initialize the segmentation of abdominal organs in [2] and [3]. Both methods used measures of relationship and hierarchy between organs and manual landmarks. Meanwhile, in [4], Zhou et al. constructed the liver atlas for improved liver segmentation that controlled the points of diaphragm surface according to TPS. Linguraru et al. [5]
suggested to segment the liver and spleen from the CT images using the position of the xiphoid to register the probabilistic atlas based on affine transformation.
Yamaguchi et al. [6] applied the region-growing method and the probabilistic at-las to segment the liver from the CT images by registered using the anatomical landmarks of the liver based on TPS. Notably, the construction of these abdominal atlases required manual landmarks through user interaction. On a different note, Okada et al. [7, 8] developed a hierarchical statistical atlas of the liver normalized to the abdominal cavity as part of a process of automatically segment the liver.
Recently, Li et al. [9] also proposed an automated liver segmentation method us-ing the rib cage to map the probabilistic atlas onto the input volume. They only considered one type of rib cage as the reference for registration.
However, there are two limitations to these segmentation methods: (1) Errors in landmark extraction may cause incorrect estimation of the transformation; (2) Even though the transformation of landmarks can be accurately estimated, it is difficult to achieve the accurate transformation of the organ because of the large variation of anatomical structure.
Taking these two points into consideration, we propose a template matching framework based on a probabilistic atlas for organ segmentation. In our proposed method, we firstly need to construct a bounding box around the organ based on human anatomical structure, which complies to the relatively rigid characteristics of organ position; and then the probabilistic atlas is used as a template to find the organ in this bounding box by the use of template matching.
Additionally, as described in [1], the main disadvantage of the conventional probabilistic atlas is that it may cause a bias for the specific patient studies due to the single reference. Numerous investigations showed that an iterative atlas construction can reduce the dependence on the reference patient [9]. Hence, we construct an iterative probabilistic atlas for our proposed method to address this bias issue.
his-Figure 4.1: The framework of our proposed segmentation method.
togram of the organ’s intensities, and a probabilistic atlas of the organ. The organ bounding box is employed to estimate the approximate location of an organ for the observed data. The histogram of the organ’s intensity is used to build the Gaussian intensity model (likelihood map). The likelihood map together with the probabilistic atlas are employed in the template matching technique. Template matching aims to estimate an initial segmentation of the organ. In the estimation of the organ bounding box, construction of the probabilistic atlas and template matching, the affine registration method is employed. Furthermore, the ’Geodesic Active Contour’ (GAC) algorithm [10] is used to refine the initial segmentation re-sults so that they can guarantee the smoothness of the final result. The framework of our proposed method is shown in Fig. 4.1.
The basic idea of our proposed method is based on Bayes framework. To estimate the organ labelLwe can follow the principle of the maximum a posteriori probability (MAP) estimation. The posteriori probability can be estimated as:
pi,j,k(c|x)∝pi,j,k(x|c)×pi,j,k(c) Liver c= 1
Spleen c= 2
(4.1) where, xi,j,k and ci,j,k are the intensity values and organ labels of voxel (i, j, k) in the CT volume, where (i, j, k) is the coordinate for each voxel. pi,j,k(x|c) is the probability of voxel (i, j, k) in the likelihood map for c-th organ, pi,j,k(c) is the probability of voxel (i, j, k) in the iterative probabilistic atlas for c-th organ, which is used as a prior probability. In the following sections, we describe how to
Figure 4.2: Schematic of organ bounding box construction.
segment the multiple organs in the abdominal region.
4.2.2 Organ Bounding Box Construction
The development of the bounding box-based method is a new trend for organ segmentation that allows easy extraction from CT volumes. An abundance of literature has been published on the bounding box for the estimation region of organs [11, 12, 13, 14, 15, 16, 17]. These methods can automatically detect the localization of organs within abdominal CT scans and can be summarized as fol-lows: regression-based [11, 12, 13], classification-based [14, 15], registration-based [16, 17]. The first two methods usually need to deal with large amounts of data and consequently need the longest computation times for the training process.
However, the registration-based approaches are more efficient and can rough-ly align the organs in onrough-ly a few seconds. For registration-based approaches, atlas-based methods have enjoyed much popularity for their conceptual simplicity.
Our algorithm incorporates the atlas concept within a bounding box construction.
Since the bone of the body is considered as a rigid part, we perform affine regis-tration to transform the bones in our work. A scheme of the organ bounding box construction is shown in Fig. 4.2, which contains the bone segmentation, the bone registration and the bounding box construction for one organ.
A. Bone segmentation
Several methods have been proposed for bone extraction [18, 19]. The main ap-proach of these methods was to estimate the intensity range of the bone and apply the region-growing technique to find the bone. Following a similar approach to these methods, our method’s first step is to find the abdominal region using Ot-su thresholding [20]. Then, we model the intensity range of bones ([Tlow, Thigh]) and was automatically selected a seed point. Finally, after the estimation of the
Figure 4.3: The flowchart of the bone segmentation method.
intensity range of the bones and the seed point, the region-growing algorithm [21]
is employed to segment bone and the morphological operators used to refine the results. The flowchart of our bone extraction method is shown in Fig. 4.3.
In order to reduce the computation demands of our algorithm, we firstly extract the abdominal region as a ROI using Otsu thresholding technique. The high threshold value (Thigh) of the bone’s intensity is defined as the largest intensity value in the abdominal region, as shown in Fig. 4.4.
We use one slice image to find the low threshold (Tlow)for bone segmentation.
Though the selection of the slice is not very critical, we select the third slice for this task because the first several slices usually do not include the kidney, which has an overlapped intensity range with the bone. In the abdominal region, there are several tissues which make estimation of the bone intensity range difficult.
However, regions that include less tissues cause estimation of bone intensities more accurate. Thus, the estimation of the bone intensity range will be more accurate.
One example of the third slice image is shown in Fig. 4.6(a). Subsequently, an iterative search approach is applied to find the bone’s low threshold, as shown in the left-hand side of Fig. 4.3. First, we assume a intensity distribution of the selected slice as a Gaussian Model (Fig. 4.4). Then we can estimate the parameters of the corresponding Gaussian distribution (G(µ, σ)). In order to segment the bone in CT data, we need to calculate a threshold value Tlow by using the parameters of the assumed Gaussian distribution. With the initial threshold value Tlow =µ
Figure 4.4: Intensity distribution of one CT volume (blue line). The yellow line is the assumed Gaussian distribution (G(µ, σ)).
(the mean value µ), we iteratively increase Tlow with a fixed value (0.2∗σ) until a ratio of the segmented bone volume to whole volume is less than a predefined value (P
). Fig. 4.5 shows the bone segmentation results with different thresholds (percentages). Lowering the threshold value for the bone segmentation results in the loss of bone regions; Conversely, with the higher threshold value for the bone, it may include multiple organs. According to our experimental results (the mean Dice measure), we set the predefined valueP
= 2% as the optimum threshold value to prevent the inclusion of the intestine, colon and blood vessel in the segmentation results in our work. By this dynamic threshold method, the low thresholdTlow for the bone is obtained.
Next, a seed point is automatically selected after the input image is thresh-olded ([Tlow, Thigh]). The whole procedure of the seed point estimation is shown in Fig. 4.6. It consists: (1) Definition of the central zone (Fig. 4.6(c)) in the ex-tracted abdominal region. We define the central zone, whose area is 25% of the abdominal region, as a ROI; (2) Generating the thresholded image (Fig. 4.6(d)) by thresholding the third slice of CT volume with [Tlow, Thigh]; (3) Obtaining the largest segmented object in the central zone; (4) Defining the central point of this largest object as a seed point (Fig. 4.6(e)). If the central point (the estimated seed point) is located within empty space, the nearest point to the central point of this largest object is defined as the seed point (Fig. 4.6(f)).
After estimating a seed point and the bone’s intensity range [Tlow, Thigh], we apply a seed region-growing algorithm [8] to segment the bone.
Figure 4.5: (a) Bone segmentation accuracy with different thresholds (P ); (b) P= 1%; (c) P
= 2%; (d) P
= 5%.
Figure 4.6: Steps of automatic extraction of a seed point in the third slice of CT volume. (a) One CT slice; (b) The extracted abdominal region; (c) The central zone of the abdomen; (d) The input image is thresholded using the optimum threshold [Tlow, Thigh]; (e) The central point (Green point) of the largest object;
(f) The final seed point.
Figure 4.7: Four different reference bones. (a) Reference 1 (shoulder to legs); (b) Reference 2 (Liver Region to legs); (c) Reference 3 (Liver Region to pelvis); (d) Reference 4 (Liver Region).
B. Bone registration
Registration is an important step in the whole process due to its primarily in-fluence on the overall accuracy of the bounding box. Since bones are considered as the rigid frameworks in the human body, the affine transform is used for the registration in our work. Empirically, the selected reference data have to approx-imately resemble the whole set of input images. However, due to variations in the abdominal region, which is scanned during the imaging procedure, the bone shape/size shows variation. Hence, only using one type of bone is not sufficient for registration. Taking this into consideration, we utilize four different reference bones for the registration. We select one training data which includes shoulder to leg region. Then we extract four different regions of the body to serve as four references, as shown in Fig. 4.7. The first reference data includes shoulder to legs.
The second reference data includes from the liver region towards legs. The third reference data includes the liver region and pelvis. The fourth data only includes the liver region.
The root mean squares (RMS) metric is used to adaptively select the appro-priate bone reference when we want to register a new data. Initially, the extracted bone for the new data is registered to all four reference bones, separately. Then, we calculate the RMS value between the new data and each reference as a measure-ment of the accuracy of the registration result. The four RMS values are compared and the extracted bone is registered to the reference bone corresponding to the lowest RMS value. Wherein, the RMS measurement is calculated based on the binary image of the segmented bone.
Figure 4.8: The bounding box for the training data. (a) Examples of bone segmen-tation; (b) Registration of bone (Gray-reference bone, Blue-moving bone before registration, Green-the moving bone after registration); (c) The bounding box of the liver (the ensemble of six registered livers); (d) The bounding box of the spleen (the ensemble of six registered spleens)
C. Organ bounding box construction
Considering the variations in anatomical structure, the positions of organs will differ for all the data. Hence, to build a bounding box for an organ, we use the various organ data as training data. Firstly, we extract the bone from one data.
Next, using the adaptive selection of the appropriate reference bone, the extracted bone is registered to the corresponding reference by using an affine transform.
Then, the same transformation is used to transform the segmented organ mask of this data. Finally, after registration of all the training data, we find the volume of interest (VOI) of the ensemble of the organs and use it as the organ bounding box.
Fig. 4.8 illustrates an example of the spleen bounding box based on six samples.
For the test data, according to the adaptive selection of bone reference, once the reference bone is registered to the bone of the observed data, the organ bounding box is automatically transformed onto this observed data using the same trans-formation. Thus, the candidate organ region is extracted based on this registered organ bounding box.
Figure 4.9: Building the spleen (L= 2) likelihood map.
4.2.3 Gaussian Intensity Model
A Gaussian intensity model (or likelihood map) of an organ refers to an image in which the intensity of a voxel corresponds to the probability of belonging to that organ. To reduce the irrelevant neighboring tissues from the input candidate organ region and to improve the segmentation accuracy, the Gaussian model is applied for analyzing the input organ intensity distribution, and can be described as:
pi,j,k(x|c) = 1
p2πσc2exp(−(xi,j,k−µc)2
2σc2 ) (4.2)
The Gaussian parameters mean µc and variance σc2 can be estimated by ob-serving the density histogram of the organ. In this work, c= 1 is the liver and c= 2 is the spleen. We only construct the gaussian intensity model for the organ of interested. We do not model the background tissues since it includes several tissues which cannot be modeled by a single Gaussian distribution. A Gaussian distributionN(µ, σ) which has the same FWHM (full width half maximum) with the observed density histogram is decided [4]. In dealing with the large range of organ intensities, we follow a non-parametric scheme and use training data to build the histogram of the voxels’ intensities. For the gray scale distribution of an organ, the 100 bin histogram was made to closely approximate the true distribution.
We convert the intensity value (xi,j,k) of voxel (i, j, k) into a probability (a likelihood probability pi,j,k(x|c)) for the organ (c) using Equation 4.2. A typical case of the spleen is shown in Fig. 4.9, in which the original image is converted into a likelihood map for spleen. Comparing the original CT image (Fig. 4.9(a)) with the likelihood map (Fig. 4.9(b)), reveals that the spleen can be more easily distinguished from other tissues in the likelihood map, since the spleen region has a higher probability. Additionally, it can be found that distinguishing the spleen from several tissues with similar intensity is difficult by depending on the intensity information (likelihood).
Figure 4.10: The diagram of the iterative atlas construction.
4.2.4 Iterative Probabilistic Atlas Model
An organ probabilistic atlas can serve as a prior probability image of the organ.
The segmented organ masks in the training data are used to construct a probabilis-tic atlas for the candidate organ. In conventional construction of a probabilisprobabilis-tic atlas, a single training data is arbitrarily chosen as the reference image, and the rest of the data are considered as the moving images and registered to this ref-erence. However, as identified by Park et al [1], the main disadvantage of the conventional atlas is that a single reference is limited in representing the whole population of potential target cases. Hence, only mapping to a single reference may cause a bias toward that specific patient study. According to surveys on this problem, this issue can be addressed by an iterative atlas construction that reduces the dependence on the reference patient [1, 9].
The construction process of our iterative atlas is illustrated in Fig. 4.10. Reg-istration of the organ masks is the most important step in the atlas construction process, which is performed using a rigid scheme in our method. In the initial step we register the organ mask in different CT cases onto the reference volume, and then the resulting atlas from the previous iteration is used as the reference for the following construction phase. The mapping result is the atlas that shows the probability of a given voxel belonging to the organ. Hence, the atlas is a compos-ite and reflects variations in organ shape and size. To accurately build the atlas, although we do not use any bone as the reference image, the corresponding bone still needs to be used to normalize the organ mask in the same physical coordinate
since the spacing and origin of training data are not the same.
In our proposed segmentation method, the probabilistic atlas is used as a tem-plate to detect the organ in the organ bounding box by the use of temtem-plate match-ing technology. Hence, when the iterative probabilistic atlas is constructed, we select the largest region (VOI) containing all the organ masks as our atlas tem-plate p(c).
4.2.5 Template Matching Framework
Template matching is one of the key innovations of our algorithm. Template matching is based on the idea that those voxels which belong to the organ have similar intensities to the organ’s histogram and their locations are the parts of the image where the organ may be found with a higher probability. The aim of this step is to build a probability image, which assigns to each voxel the probability of belonging to the organ of interest. To build such a probability image, we use two probability images: a gaussian intensity model (a likelihood map) and a probabilis-tic atlas. Our algorithm incorporates a template matching concept within a Bayes framework, according to the definition of Equation 4.1. p(x|c) is the likelihood map which can be obtained using the organ histogram and p(c) is estimated by using the organ probabilistic atlas. We use the constructed organ bounding box, described in previous step, to restrict the processing region in template match-ing. The size of likelihood map is the same as the original CT image. However, considering the probabilistic atlas can be regarded as a template in the matching process, the size of probabilistic atlas is usually smaller than the likelihood map.
The search range for atlas matching is the bounding box, in which the organ is located. A sub-volume that has the same size and location with the probabilistic atlas is extracted from the likelihood map. Then a probability image is generated by voxel-wise multiplication of the sub-volume with the probabilistic atlas. For this generated probability image, the total probabilities (the summation of the probability image) are used as a similarity measure. We search from the bottom-left voxel (start point) of the bounding box to the top-right voxel (end point) of the bounding box to find the best object location (dopt, eopt, fopt) that maximizes the similarity (the summation of the probability image). At this location the sub-volume of the likelihood map will reach the highest matching with the probabilistic atlas. This procedure can be formulated as:
(dopt, eopt, fopt) = arg max
(d,e,f)∈ΩB
"
X
i,j,k
pi+d,j+e,k+f(x|c)×pi,j,k(c)
#
(4.3) where ΩB is the bounding box for the organ and (d, e, f) represents the location of atlas in the bounding box ((a, b, c)∈ΩB). Furthermore, the generated
probabil-ity image at location (dopt, eopt, fopt) is padded into the corresponding location of a empty volume that has the same size with the likelihood map. This volume serves as the output of template matching, which contains trivial non-organ structures.
In order to remove trivial objects, we firstly threshold this output volume and then smoothed its surface using an opening morphology operator. The generated results is the atlas-based segmented organ, which can be regarded as the initial segmentation results.
4.2.6 Refinement of Organ Segmentation
The final step of organ segmentation is to refine the initial segmentation result by using a Geodesic active contour (GAC) algorithm [10].
GAC is an active contour method which tries to solve a differential equation using an iterative method. The differential equation is the result of minimiza-tion of an energy funcminimiza-tion which moves a contour toward the true boundary of the object. For the GAC segmentation algorithm, an initial contour is required.
Its segmentation results are sensitive to the initial position. As explained in the previous section, after template matching, we can achieve an atlas-based segment-ed organ. Hence, in our work, the generatsegment-ed atlas-bassegment-ed segmentsegment-ed organ was used as an initial contour of the GAC algorithm to achieve a more accurate final segmentation result.