Segmentation of cytoplasm and nuclei of abnormal cells in cervical cytology using global and local graph cuts.

G Model

ARTICLE IN PRESS

CMIG-1241; No. of Pages 12

Computerized Medical Imaging and Graphics xxx (2014) xxx–xxx

Contents lists available at ScienceDirect

Computerized Medical Imaging and Graphics journal homepage: www.elsevier.com/locate/compmedimag

Segmentation of cytoplasm and nuclei of abnormal cells in cervical cytology using global and local graph cuts Ling Zhang a,b,c,d , Hui Kong e , Chien Ting Chin a,c,d , Shaoxiong Liu f , Zhi Chen b , Tianfu Wang a,c,d,∗ , Siping Chen a,c,d,∗ a

Department of Biomedical Engineering, Shenzhen University, Shenzhen 518060, China Department of Electrical and Computer Engineering and Iowa Institute for Biomedical Imaging, The University of Iowa, Iowa City 52242, USA c National-Regional Key Technology Engineering Laboratory for Medical Ultrasound, Shenzhen 518060, China d Guangdong Key Laboratory of Biomedical Information Detection and Ultrasound Imaging, Shenzhen 518060, China e Department of Mechanical Engineering, Massachusetts Institute of Technology, Cambridge, MA 02139, USA f Department of Pathology, People’s Hospital of Nanshan District, Shenzhen 518052, China b

a r t i c l e

a b s t r a c t

i n f o

Article history: Received 6 August 2013 Received in revised form 3 February 2014 Accepted 11 February 2014 Keywords: Abnormal cells Cervical cell segmentation Graph cut-based segmentation Touching-nuclei splitting

Automation-assisted reading (AAR) techniques have the potential to reduce errors and increase productivity in cervical cancer screening. The sensitivity of AAR relies heavily on automated segmentation of abnormal cervical cells, which is handled poorly by current segmentation algorithms. In this paper, a global and local scheme based on graph cut approach is proposed to segment cervical cells in images with a mix of healthy and abnormal cells. For cytoplasm segmentation, the multi-way graph cut is performed globally on the a* channel enhanced image, which can be effective when the image histogram presents a non-bimodal distribution. For segmentation of nuclei, especially when they are abnormal, we propose to use graph cut adaptively and locally, which allows the combination of intensity, texture, boundary and region information. Two concave points-based approaches are integrated to split the touching-nuclei. As part of an ongoing clinical trial, preliminary validation results obtained from 21 cervical cell images with non-ideal imaging condition and pathology show that our segmentation method achieved 93% accuracy for cytoplasm, and 88.4% F-measure for abnormal nuclei, outperforming state of the art methods in terms of accuracy. Our method has the potential to improve the sensitivity of AAR in screening for cervical cancer. © 2014 Elsevier Ltd. All rights reserved.

1. Introduction Screening by cytology is the commonest approach to prevent cervical cancer at a pre-cancerous stage [1]. However, it is known that screening of cervical cytology slides is labor intensive and mentally demanding on the cytotechnologist [2]. Automation-assisted reading (AAR) techniques have the potential to increase productivity and reduce screening errors. AAR uses automated microscopy to collect images of cervical cells, and then performs cell segmentation to quantitate cellular structures and morphology in order to extract candidate cells for targeted reading by pathologists. A large, prospective randomized trial found that although AAR increased the productivity by 60–80%, the sensitivity was reduced [3]. To

∗ Corresponding authors at: Department of Biomedical Engineering, Shenzhen University, Shenzhen 518060, China. Tel.: +86 755 86671932. E-mail addresses: [email protected] (T. Wang), [email protected] (S. Chen).

improve the sensitivity, it is important to segment the abnormal cells reliably. A variety of segmentation methods for cervical cytology have been proposed in recent years. The majority of cytoplasm segmentation used one or multiple of the following techniques: K-means [4,5], edge detection [6], thresholding [7,8] and active contours [5,9]. Most of these works are designed for images of isolated cells, especially for those in the Herlev data set [11]. For segmentation in images containing multiple cells, thresholding [7,8,12] and level set [9,10] techniques have been used. However, thresholding may lead to inaccuracy due to non-ideal imaging conditions, such as inconsistent staining and/or illumination and overlapping cells. The level set method is also suboptimal, due to computation cost and its tendency to locate some local extrema (e.g., inhomogeneous illumination). For nucleus segmentation, related works can be divided into three groups: single-nucleus segmentation, multiple-nuclei segmentation, and touching-nuclei splitting. Existing singlenucleus segmentation approaches utilize contour, shape and color

http://dx.doi.org/10.1016/j.compmedimag.2014.02.001 0895-6111/© 2014 Elsevier Ltd. All rights reserved.

Please cite this article in press as: Zhang L, et al. Segmentation of cytoplasm and nuclei of abnormal cells in cervical cytology using global and local graph cuts. Comput Med Imaging Graph (2014), http://dx.doi.org/10.1016/j.compmedimag.2014.02.001

G Model CMIG-1241; No. of Pages 12

ARTICLE IN PRESS L. Zhang et al. / Computerized Medical Imaging and Graphics xxx (2014) xxx–xxx

2

information by using active contour model [5,13], parametric fitting [14] and difference maximization [4]. Nucleus segmentation for images containing multiple cells employ thresholding [9,15], Hough transform [7], morphology (watershed) [8,16,17], and level set [10] techniques. Some of these methods [7–9] assume all nuclei located within the cytoplasm, as a result in that bare nuclei, which are strongly indicative of early cervical cancer, are missed. For splitting of touching-nuclei, morphological erosion [15] is a commonly used approach. An unsupervised method [18] for splitting touching nuclei is proposed by a combination of distance transform, expectation–maximization algorithm and ellipse fitting techniques. Through representing the nuclear shape by the vibrations of a spring-mass system, and learning the vibration models by active shape model, the nuclear boundary in the overlapping areas can be obtained accurately [19]. Generally, there are many cervical cells with mutual overlaps in a field-of-view (FOV) (see Fig. 1(a) as an example). Hence, a segmentation algorithm should capture the cytoplasm, multiple nuclei and touching-nuclei. Three of the aforementioned methods [7,8,10] can achieve all of these three tasks. However, these methods were developed to segment images, which only contain healthy cells rather than the mixture of healthy and abnormal cells. Since the size, shape and chromatin distributions of abnormal nuclei vary significantly (Fig. 1(b)), further development is needed to address typical variations encountered in a clinical setting. So far we have found only one previous study of automated segmentation of abnormal nuclei [20] and only very limited data was reported. Segmentation of nuclei in histological images employ adaptive thresholding combined with active contour models [21]. Such automated methods are comparable with the manual delineation in segmentation accuracy. Similarly, graph cut (GC) approaches [22] are highly attractive in nucleus segmentation. The binarization of nuclei based on GC is addressed in [23] with results being more accurate than global thresholding. Prior knowledge like nucleus shape [24], manual annotation and local image features [25] can be incorporated in the GC framework to allow more robust segmentation. In this paper, we propose to segment cervical cytoplasm globally and segment cervical nuclei locally based on GC approaches, respectively. The overall method, the outline of which is shown in Fig. 2, follows two tracks. The first track applies cytoplasm segmentation method to separate the cytoplasm from the background (upper part in Fig. 2). The second track applies nucleus binarization and touching-nuclei splitting methods to segment nuclei (lower part in Fig. 2). Specific contributions of the presented work consist of: (1) cytoplasm segmentation by using the multi-way GC [22] on the a* channel (CIE L*a*b* model in Ref. [26]) enhanced images (Section 2.1); (2) nucleus segmentation using a proposed “Local Adaptive GC” (LAGC) method which results in robust binarization of nuclei with a variety of morphologies, chromatin distributions and low contrasts (Section 2.2); (3) a splitting method for touchingnuclei by combining two concave points-based methods [27,28] (Section 2.3).

2. Methods 2.1. Cytoplasm segmentation Given a cervical cell image, the a* channel in CIE LAB color space is used for preprocessing. Initial segments are generated automatically by using Otsu’s multiple thresholding algorithm [29] on the preprocessed image. The segmentation is refined by a multi-way GC method. In the rest of this section, we introduce the details of our cytoplasm segmentation.

2.1.1. Pre-processing In cervical cell images, the poor contrast, non-uniform staining and noise will likely hinder cell segmentation. To enhance the contrast, we extract the a* channel in the CIE L*a*b* color space [26]. The reason for selecting the a* channel is as follows: In Hematoxylin and Eosin (H&E) staining, cell regions are colored with tones of red and background regions remain colorless. This inspires us to use color to discriminate cells from background. The a* channel represents change between red and green and is able to embody this difference. Hence, in the a* channel image, cells are obviously brighter than the background. To further enhance contrast, the a* channel image is stretched linearly from their original intensity range [Imin , Imax ] to the range [0,255]. A common and effective technique in handling noise is the median filter. It was demonstrated by Tsai et al. [4] that median filter can eliminate both impulse and Gaussian noise in cervical smear images. In our work, a 5 × 5 median filter is applied to the contrast enhanced images to discard noise. Fig. 3(b) shows the result of the above preprocessing.

2.1.2. Global segmentation (multi-way GC) Although the difference between cells and background is significantly enhanced by preprocessing, not all the image histograms present bimodal distribution due to the complexity of our images which include inhomogeneous illumination, non-uniform staining and the presence of inflammation cells and debris. Therefore, a single threshold cannot successfully separate the cervical cell from the background. For example, some cytoplasm with brighter intensity tends to be classified into the background class. With this in mind, we propose to use a multi-way GC approach. The output image contains four classes. The class with the lowest mean intensity is the background. The other three classes, which contain cytoplasm, nuclei, inflammation cells and debris are integrated as the foreground (cell region). Our approach is summarized below: • Given an image with intensity values from Imin to Imax , we compute three optimal thresholds th1 *, th2 * and th3 * with Imin ≤ th1 * ≤ th2 * ≤ th3 * ≤ Imax , by applying 3-thresholds Otsu algorithm [29]. Then, the mean intensity values, c1 , c2 , c3 and c4 , of the four classes (C1 , C2 , C3 and C4 ) are computed, where C1 = [Imin , . . ., th1 *], C2 = [th1 * + 1, . . ., th2 *], C3 = [th2 * + 1, . . ., th3 *], and C4 = [th3 * + 1, . . ., Imax ]. • With c1 , c2 , c3 and c4 , we construct a four-terminal graph, and employ the Potts model [22] energy function Ep (f) EP (f ) =

p∈P

Dp (fp ) +

ω{p,q} · T (fp = / fq ),

(1)

{p,q} ∈ N

where f denotes the pixel label, p indexes pixels, N is the set of adjacent pixels, {p, q} represents a pair of pixels. Potts model is a metric: V(˛, ˇ) = K·T(˛ = / ˇ), where T(·) is 1 if its argument is true, and 0 otherwise. In this case, discontinuities between any pair of labels are penalized equally. The Potts model parameter K is used to control the tradeoff of weights between the data term and the pixel continuity term. Higher value of K strengthens the neighborhood constraints and makes the segments larger, and vice versa. Boykov et al. [22] demonstrated that the minimization of Ep (f) can resolve the multi-way cut problem. The first term of Eq. (1) is the data term, which is determined by the connection energy t-links between each pixel and each terminal of the graph. In our work, Dp (fp ) is assigned as (cp − Ip )2 , where Ip is the intensity value of pixel p. The second term in Eq. (1) is the pixel continuity term, which is determined by the connection energy n-links between


G Model

ARTICLE IN PRESS


L. Zhang et al. / Computerized Medical Imaging and Graphics xxx (2014) xxx–xxx

3

Fig. 1. (a) A sample of cervical cell image with boundaries of the cytoplasm and the nuclei marked as yellow and green, respectively. Some abnormal, healthy and overlapped nuclei are marked by rectangles with red, light blue and black colors, respectively. (b) Zoomed region according to the red rectangles in (a). (c) Zoomed region according to the light blue rectangles in (a). (d) Zoomed region according to the black rectangles in (a).

neighborhood pixels in the multiple-terminal graph. The ω{p,q} is defined as in [22]:

ω{p,q} = ω(Ip − Iq ) =

2K

if |Ip − Iq | ≤ 5

K

if |Ip − Iq | > 5

function of Eq. (1) is optimized, and the global optimal pixel category label f is obtained. Consequently, we can re-label each pixel to its new class.

(2)

In our work, the Potts model parameter K is set as 0.1 through empirical evaluation. Ip and Iq are two neighboring pixel values. • Finally, by using an implementation of the ␣-expansion and fast max-flow/min-cut algorithms introduced in [22,30], the energy

Finally, to split some cytoplasm with slight overlaps, the segmented binary image (after merging the three foreground subclasses) is processed by morphological opening with a 5 × 5 square structuring element. Actually, there exists some other method to replace the Otsu’s multiple thresholding, such as K-means and Gaussian Mixture Models (GMM) clustering [31]. We choose Otsu


G Model CMIG-1241; No. of Pages 12 4


Fig. 2. Block diagram of overall approach.

thresholding because it is more efficient and shows better segmentation results in our evaluation (we here refrain from these results). In addition, Otsu’s method is one of the best thresholding methods for general images [32]. Fig. 3(c)–(e) shows the examples of cytoplasm segmentation. By examining the center part of Fig. 3(f), we see that the segmentation is not affected by the bright illumination and dirt.

2.2. Nucleus segmentation Since the presence of bare nuclei is a strong indicator of disease, our nucleus segmentation does not depend on cytoplasm segmentation. The nucleus segmentation consists of a nucleus binarization algorithm and a touching-nuclei splitting method. For the nucleus binarization, we use GC approach adaptively and locally, denoted

Fig. 3. (a) An original color image. (b) Contrast enhancement and noise removal. (c) 3D Otsu. (d) Multi-way GC. (e) Foreground merging and morphological open. (f) The obtained cytoplasm boundary.


G Model

ARTICLE IN PRESS



5

as LAGC. LAGC first uses an efficient adaptive thresholding [33,34] to detect the nuclei and then refine each nucleus within its local neighborhood using GC [22,30]. For the touching-nuclei splitting, given a connected component of the binary component, we initially estimate whether it is a touching-nuclei clump or not based on morphological and gradient features, and then combine two concave point-based algorithms [27,28] for splitting.

corresponding to a larger rectangle (with size [slength + s , swidth + s ]) is first extracted from the preprocessed image [Fig. 4(b) and (c)]. Then the subimage is stretched linearly to further enhance the contrast. The histograms of subimages are found to be bimodal and can be modeled well by a mixture of two Poisson distributions. This modeling choice is based on the analysis in Ref. [23]. The Poisson mixture parameters are given:

2.2.1. Pre-processing Before applying the nucleus segmentation algorithm, a practical problem should be considered: if the cytoplasm is deep-stained or some inflammation cells cluster together, there might be many abnormal segments, which are often the major challenge in cervical cytology automation [35]. To enhance the contrast between the nuclei and cytoplasm, the original color image is preprocessed using the procedure designed in our previous work [28]. Briefly, we convert the original RGB color space to HSV, in which the V channel is extracted and enhanced through linear stretching. Finally, the median filter with a mask of size 5 × 5 pixels is applied to discard noise. After the preprocessing, the nuclei tend to become much darker while the cytoplasm is much brighter.

0 =

2.2.2. Local segmentation (nucleus binarization) Considering the sensitivity requirement of an AAR system, we propose a GC based nucleus binarization method which works in an adaptive and local paradigm, as illustrated in Fig. 4. The adaptive stage detects each nucleus region approximately [Fig. 4(a) and (b)] by applying an efficient adaptive thresholding algorithm that uses intensity and texture information. The local stage refines each adaptive segment [Fig. 4(c) and (e)] within its local neighborhood by using a Poisson distribution based GC which utilizes boundary and region information. This approach can not only overcome the influence of non-uniform staining and illumination, but can also improve the binarization of nuclei with non-uniform chromatin distribution and with low intensity difference vs. the surrounding cytoplasm. • Adaptive stage. Because of the inhomogeneous illumination and non-uniform staining, it is hard to define a global threshold without either missing some of the nuclei or segmenting non-nuclei parts. Therefore, adaptive thresholding is desired. Trier and Jain [36] evaluated 11 thresholding algorithms on document images, and concluded that the Niblack algorithm performed the best. Niblack algorithm has also found applications in the segmentation of cell nuclei [21]. Sauvola et al. [33] proposed a new adaptive binarization method, and the benchmarking results showed that their method outperformed the others including the Niblack algorithm. Therefore, in our work, we use this method to compute a threshold t(x, y) for each pixel at location (x, y):

t(x, y) = m(x, y) 1 + k

s(x, y) R

−1

,

(3)

where m(x, y) and s(x, y) are the mean (intensity information) and the standard deviation (texture information) of the gray level values within a w × w pixel window of (x, y). The R is the maximal value of standard deviation and set to 128, and k is a constant which takes values in the range [0.2, 0.5]. If the intensity value of a pixel is lower than t(x, y), it is segmented as a nucleus pixel. Parameter w can be chosen in relation to the maximal size of the nuclei (71 in our work). A k-value of 0.3 gives good results in our work. With the integral-valued images method [37], this method runs very fast [34]. • Local stage. This stage aims at achieving more accurate binarization of the nuclei, especially the abnormal nuclei with non-uniform and weak staining. Based on the bounding box of each adaptive segment (with size [slength , swidth ]), a subimage

1 Nfg

Ip × h(Ip ),

1 =

Ip ∈ fg

1 Nbkg

Ip × h(Ip ),

(4)

Ip ∈ bkg

where 0 and 1 are the mean intensity values of the foreground and the background, respectively. Nfg and Nbkg are the number of pixels in the foreground and background, respectively. h(Ip ) represents the normalized image histogram, NIp

h(Ip ) =

Nfg + Nbkg

,

(5)

where NIp is the number of pixels with intensity value Ip . In some subimages, a weakly stained nucleus may be surrounded by some dark objects (e.g., nuclei or artifacts, note the gray areas in the first column of Fig. 4(c) as examples). Taking these dark objects as the foreground may lead to an inaccurate segmentation result. Hence, in each subimage (with labeled foreground), we exclude objects which neither belong to the foreground nor the background. Then, 1 of each subimage is used to replace the intensity values of the excluded objects. With 0 and 1 , we compute Poisson probabilities of the foreground and background as follows: fg(Ip ) =

u0 Ip e−0 , Ip !

bkg(Ip ) =

u1 Ip e−1 , Ip !

(6)

With these probabilities, we construct the following GC energy function: E(f ) =

p∈P

Dp (fp ) +

Vp,q (fp , fq ),

(7)

{p,q} ∈ N

where the first term has two possible values depending on whether the foreground or background model is used, and the second term is the pixel continuity term. They are written as follows, Dfg (fp ) = − ln fg(Ip ), Dbkg (fp ) = − ln bkg(Ip ),

Vp,q (fp , fq ) = ı(fp , fq ) · exp

−

(Ip − Iq ) 2 2

2

,

(8) (9)

where ı(·) is 0 when fp = fq , and 2 otherwise. is the scale factor which is determined by the smoothness of nuclear chromatin, and can be set to [2 0 4 0]. In our work, it is set to 30. The above energy function is represented as a weighted graph through GC. By seeking the minimum cut of the graph through the maxflow/min-cut algorithm [22], the optimal segmentation of the image is obtained. Note that Al-Kofahi et al. [23] also used Poisson distributions in their GC-based binarization. The difference is that the Poisson distributions are used to seek a global optimal threshold to separate the foreground and the background for the GC initialization in [23], while the foreground and background are already known after the adaptive stage in our method. 2.2.3. Post-processing After the refinement process (local stage using GC), for all segments in the subimage, only the segment which has the maximal overlap area with the object in the adaptive stage is retained. To reduce the computational complexity, for a certain adaptive segment (Eq. (2)), we empirically set a condition based on one feature, namely, roundness Fr . The Fr , proposed by Haralick [38],


G Model

ARTICLE IN PRESS



6

Fig. 4. (a) The input image. (b) The result after preprocessing and adaptive thresholding. The white objects need the LAGC process while the gray objects needn’t. (c) Refine some of the adaptive segments within its local neighborhood by using a Poisson distribution based GC. (d) Replace the adaptive segments by the refined results. (e) LAGC binarization result.

is defined as: Fr = R / R , where R and R are the means and variance of the distance between the centroid and the boundary points of , respectively. Then, the which satisfies Fr < 1.2 or Fr > 3 should be the well segmented nuclei or the artifact cluster, but not wrongly segmented nuclei, and it needs no further refinement. The gray objects in Fig. 4(b) are some examples of the adaptive segments which need no refinement. 2.3. Touching-nuclei splitting and reconstruction The ability of splitting touching nuclei is crucial to a fully automated cervical nucleus segmentation method. In our work, we combine two of our previously proposed concave points-based methods [27,28] to split the touching-nuclei. This combination integrates morphological features (geometric center and arc chord ratio) and gradient feature (radial symmetry center [39]). Briefly, a certain connected region i is deemed as a touching-nuclei clump if it satisfies the following two conditions: (Fs > 1.1 and Fr < 4) or |ri − gi |2 > 5,

(10)

where ri is the most likely radial symmetry center, gi is the geometrical center, and |.|2 is the Euclidean distance. Fs is the shape factor and is defined as: Fs = L2 /4␲Fa , where L is the perimeter and Fa is the area. If i only satisfies the former condition, we use the arc chord ratio based algorithm [28] to split i instead of empirically reallocating the geometrical center as in [27]. Otherwise, i is split by using the radial symmetry based splitting algorithm in [27]. The splitting procedure is illustrated in Fig. 5(b) and (c). Note that some debris and inflammatory cell clumps are eliminated by thresholding their sizes and shape factors. The splitting lines obtained by the above method usually cannot accurately delineate the occluded contour, and thus will influence the reliability of nuclei feature extraction. In our work, the constrained ellipse fitting technique [18] is employed for occluded contour reconstruction. Fig. 5(d) shows some reconstructed occluded contours. 3. Experimental methods 3.1. Clinical data collection All images used in this study were acquired using an Olympus BX41 microscope equipped with 20× objective (Olympus America,

Inc., Central Valley, PA), Jenoptik ProgRes CF Color 1.4 Megapixel Camera (Jenoptik Optical Systems Inc., Jena, Germany), and MS300 motorized stage (NJRGB Inc., Nanjing, China). Image specifications were 24 bit RGB channels with resolution of 1360 × 1024 pixels. The data set included 51 cervical cell images from 21 cervical slides, which were collected from the Department of Pathology, Shenzhen Sixth People’s Hospital Huazhong University of Science and Technology Union Shenzhen Hospital, Shenzhen, China, in 2010. All slides are prepared using manual liquid-based cytology (MLBC) technique [40], and are stained with H&E. Out of 51 images, 30 were used to tune the parameters of the algorithm and the validation was performed on the remaining 21 images. Among the 21 test images, 15 contained abnormal cells and the other 6 were normal. The ground truth used for evaluations of cytoplasm and nucleus segmentation was obtained by manual delineation by a pathologist. Only a subset of nuclei which can be unambiguously determined by human expert was annotated. The non-annotated nuclei have almost no effect on the evaluation result, since the ambiguously determined nuclei will not be checked in clinical screening by pathologists. 3.2. Quantitative assessment methodologies The evaluation of cytoplasm segmentation was based on comparison with three other cytoplasm segmentation methods [8,9,12] using Dice similarity coefficient (DSC) [41], which is defined as follow, DSC =

2|RGT ∩ RSeg | , |RGT | + |RSeg |

(11)

where RGT denotes the ground truth region, RSeg the segmented region, and |·| the number of pixels in a certain region. The results of our nucleus binarization method were compared with the output of mean shift (MS) [42], Active contours without edges (ACWE) [43], Li et al.’s method [5], and Al-Kofahi et al.’s method [23] using a pixel-based and an object-based criteria. In the pixel-based criterion, the precision (prec), recall (rec) and F-measure (F) were used as performance indices as in [7,8,10,25]. The prec and rec indicate the fraction of the amount of nucleus correctly identified in the segmented object, and in the reference ground truth, respectively; F measures the harmonic mean of prec and rec. In addition, a measure named overlap (overl) is used. These indices




7

Fig. 5. An example of touching-nuclei splitting result. (a) An original image patch. (b) Binarization image overlaid with radial symmetry centers (green dot), geometrical centers (blue cross), and concave points (red dot). (c) Splitting result. Red represents objects with no need for splitting. Although some green objects need splitting, the concave points cannot be found on their boundaries. (d) Final binarization, splitting and reconstruction boundaries overlaid on the original image. Note that we only reconstruct overlapping nuclei which have an area higher than 400 pixels (nucleus with area smaller than 400 pixels must be normal nucleus, so do not need accurate reconstruction in cervical cancer screening task). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

are defined as follows: TP prec = , TP + FP overl =

TP rec = , TP + FN

since it is one of the most important nuclear features to distinguish abnormal and healthy nuclei [46].

2 · prec · rec F= , prec + rec

TP , TP + FN + FP

4. Results and discussion (12)

where TP represents the pixel numbers of nuclei correctly identified, FP denoting the number of pixel in background which is wrongly identified as nuclei, and FN being the number of pixel in nuclei missed by segmentation. These indices can be used to evaluate all the nuclei in the image, but also to evaluate each abnormal nucleus individually. The pixel-based criterion serves as a basis for designing the object-based criterion. Nucleus detection was considered as a true positive if its recall was higher than 60%. In addition, to further estimate the binarization accuracy of abnormal nuclei, we graded the F-measure values in ranges of 0–75% (poor), 75–100% (acceptable), and 90–100% (very accurate), and then calculated the percentage of the number of abnormal nuclei in each range. To evaluate the nuclei splitting accuracy, three types of errors were employed in accordance with [23,27]: the under-splitting error, over-splitting error, and encroachment error. Furthermore, the pixel-based precision, recall and F-measure were used to evaluate the performance of splitting and reconstruction of overlapped nuclei. 3.3. Statistical correlation analysis and reproducibility analysis Linear regression analysis [44] and Bland–Altman plots [45] were used to evaluate the relationship and agreement between the manual and automatic segmentation. Considering that the segmentation of abnormal nuclei is the most important issue to an automated screening system, we asked the pathologist to manually segment the abnormal nuclei in the test images again after more than one year. Nuclear area is used as the evaluation measure

4.1. Parameters tuning Our methods were implemented using C++ on a 64-bit Windows PC, which has a 2.66 GHz quad-core CPU and 4 GB of RAM. In our method, there are two major parameters that we tuned on training images. The cluster number in the global GC was tuned on 30 training images. We set it to 3, 4, 5, and 6, respectively, The resultant cytoplasm segmentation accuracy of DSC was 0.89, 0.93, 0.93, and 0.91, respectively. The s in the LAGC was tuned on 30 nuclei images. We varied the value of s (3, 5, 10, 15). The resultant Fmeasure was 0.80, 0.84, 0.88, and 0.88, respectively. Considering the computational burdens, we suggested using cluster number = 4 in global GC and s = 10 in LAGC. Other parameters (e.g., window size w) were tuned reasonably or were set as the same as our previous studies (parameters in Eq. (10)) [27,28]. Consequently, for a new imaging protocol, these parameters can be easily adjusted by utilizing a small set of representative data. 4.2. Evaluation of cytoplasm segmentation Fig. 6 gives two examples of our cytoplasm segmentation approach on non-ideal cervical cell images, demonstrating that our method produced correct segmentations of cytoplasm in the presence of overlap cells (Fig. 6(a)), poor contrast and non-uniform staining (Fig. 6(b)). Our algorithm can also extract cytoplasm boundaries with the presence of inflammatory cells. Table 1 shows the mean and standard deviation of the cytoplasm segmentation accuracy achieved by existing algorithms [8,9,12], which are lower than our method. In Gençtav et al.’s [8] method,


G Model

ARTICLE IN PRESS



8

Fig. 6. Examples of the proposed cytoplasm segmentation method applied to two cervical cell images. (a) Heavily overlap cells. (b) Poor contrast and non-uniform staining.

Table 1 Mean and standard deviation of DSC of cytoplasm segmentation performance of the proposed algorithm compared to existing methods in Refs. [8,9,12].

Mean Std.

[8]

[9]

[12]

Ours

0.64 0.26

0.68 0.15

0.76 0.18

0.93 0.03

Table 3 Comparison of average nucleus binarization performance using object based criterion.

Mean Std.

MS

ACWE

[5]

[23]

LAGC

0.93 0.06

0.77 0.18

0.80 0.35

0.78 0.14

0.99 0.01

the only parameter is the radius of the disk of the black top-hat algorithm. We used a radius of 210 pixels as set in [8]. For Harandi et al.’s [9] method, the segmentation results will be different if the parameter and iteration step of the active contour model have different values. Based on our empirical test, the iteration step is set as 5,000 to ensure the convergence of the contour, and is set as 0.8 in accordance with [9]. Our cytoplasm segmentation method takes 0.3 s per image on average.

Table 4 Comparison of abnormal nucleus binarization performance by grading F-measure.

4.3. Evaluation of nucleus segmentation

nucleus binarization results into poor, acceptable and very accurate by using the F-measure. 93.8% of the abnormal nuclei were binarized with good quality (F-measure > 75%) by LAGC. From these three tables, we can see that the LAGC outperforms the other four methods. It should be noted that in these evaluations, overlapping nuclei are considered as a whole object. For mean shift algorithm, we directly use the EDISON software available for free download [42]. By tuning parameters on the test images, we set the spatial bandwidth hs , the color bandwidth hr and the minimum region as 10, 6.5 and 130, respectively. Intensity stretch and median filter are used before the mean shift algorithm, since we find that these preprocessing methods can bring better segmentation performance. For the ACWE algorithm, only median filter is used for preprocessing, since intensity stretch will reduce the segmentation performance in our images. We find that the change in ACWE’s parameters does not lead to significantly different segmentation results. We choose the parameters of ACWE as follows: = 1 = 2 = = ε = h = 1, t = 0.1, and the iteration step

Fig. 7 shows the comparison of nucleus binarization results on an image patch, in which there are two abnormal nuclei (green arrows). The result of LAGC is much better than those of other algorithms in terms of the binarization accuracy, especially the abnormal nucleus binarization accuracy. For example, the upper abnormal nucleus is wrongly binarized by ACWE (Fig. 7(c)), Li et al.’s method (Fig. 7(d)), and Al-Kofahi et al.’s method (Fig. 7(e)), and the lower abnormal nucleus is wrongly binarized by all the other four algorithms. Table 2 shows the comparison of LAGC and the other four algorithms in terms of average precision, recall, F-measure and overlap of binarization for both normal and abnormal nuclei. The LAGC has a 0.873 F-measure and a 0.884 F-measure on binarization of all nuclei and abnormal nuclei, respectively. Table 3 compares LAGC and the other four algorithms in terms of nuclei detection rate. The LAGC achieves a 0.99 detection rate of nuclei. Table 4 grades the abnormal

Method

% Poor

% Acceptable

% Very accurate

[5] [23] ACWE MS LAGC

45.3 48.4 39.1 25.0 6.2

54.7 51.6 60.9 75.0 93.8

17.2 20.3 25.0 34.4 51.6

Table 2 Comparison of average nucleus binarization performance using pixel based criterion. Method

MS ACWE [5] [23] LAGC

All nuclei

Abnormal nuclei

prec

rec

F

overl

prec

rec

F

overl

0.94 0.70 0.52 0.79 0.85

0.74 0.71 0.77 0.67 0.90

0.826 0.683 0.598 0.710 0.873

0.71 0.54 0.45 0.57 0.78

0.80 0.84 0.62 0.77 0.88

0.70 0.75 0.86 0.61 0.91

0.722 0.715 0.688 0.627 0.884

0.60 0.66 0.56 0.52 0.81

Note: precision (prec), recall (rec), F-measure (F), overlap (overl).




9

Fig. 7. Comparison of nucleus binarization results with (a) ground truth, (b) mean shift [42], (c) ACWE [43], (d) Li et al.’s method [5], (e) Al-Kofahi et al.’s method [23], (f) LAGC. The edges of binary masks are overlaid on the color images.

Fig. 8. Automatic segmentation results for full cervical cell images. The upper two images contain abnormal cells, while the lower two are normal cases. The boundaries of the cytoplasm and the nuclei delineated by our methods are marked as yellow and green, respectively.

is 15. The initial contour is set as a rectangle around the image. The original version of Li et al.’s method [5] is designed for single cervical cell images. We extend it to process all nuclei in an image. Our extension includes eliminating the shape constraints of candidates and performing the Radical GVF (RGVF) snake on all candidates. Following [5], the parameters ˛, , ˇ, ı, , are set as 1, 1, 5, 0.5, 10, and 2, respectively. Al-Kofahi et al.’ method is implemented in the Farsight open source project [23]. Their method includes a nucleus binarization algorithm and a touching-nuclei splitting algorithm. In comparison, only the GC based binarization algorithm is utilized.

4.4. Evaluation of touching-nuclei splitting and reconstruction Table 5 shows the performance of overlapped nuclei splitting. 3 of the 31 over-splitting nuclei are abnormal nuclei, and there are Table 5 Evaluation of our splitting performance.

Number Percentage %

Overlap nuclei

Correct split

Under-split

Encroach. errors

549 –

508 92.5

15 2.7

26 4.7




Fig. 9. Statistical correlation analysis between automated method and manual segmentation Manual 1 (a), and between manual segmentation Manual 2 and Manual 1 (b). In each graph, the upper figure is linear regression analysis results by comparing nuclear area, and the lower figure is the Bland–Altman plot.

totally 64 abnormal nuclei in our test images. That is to say, the segmentation accuracy of abnormal nuclei will suffer a 4.7% decrease due to over-splitting. In the 508 correct split nuclei, the average precision, recall and F-measure of the reconstruction results are 0.92, 0.87, and 0.89, respectively. Fig. 8 shows the final cell segmentation results from our test image set. It can be seen that our proposed methods can accurately delineate the boundaries of cytoplasm in H&E stained images in presence of inhomogeneous illumination, inconsistent staining and dirt occlusion, and achieves promising segmentation results for nuclei/overlapping nuclei having weak staining and non-uniform chromatin distribution. The average time-cost for the whole nucleus segmentation procedure is about 1.6 s per image. 4.5. Statistical correlation analysis and reproducibility analysis Fig. 9(a) shows the linear regression analysis by comparing nucleus areas and Bland–Altman plots for the automated segmentation versus Manual 1. Fig. 9(b) presents the reproducibility (intra-observer variability) assessment of manual segmentation Manual 2 versus Manual 1. These results show that: (1) the intraobserver reproducibility has a high correlation with r = 0.984. In comparison, the automated segmentation also achieves a relatively high correlation with the Manual 1 segmentation (r = 0.938). (2) The Bland–Altman plots reveal that the 95% limits of agreement were [−77,51] and [−36,22] for the automated method versus Manual 1,

and Manual 2 versus Manual 1, respectively. The automated versus Manual 1 showed a larger bias compared to the intra-observer variability.

4.6. Limitations of the reported method Our approach has some limitations. The first limitation is that it does not attempt to delineate the cytoplasm boundary for each cell. This is because reliable delineation of the cytoplasm boundary for each cell is unrealistic even for human experts in the presence of heavily overlapping cells. Furthermore, our method may fail to delineate the cytoplasm boundary in red regions caused by intensively bright illumination. Fig. 10 shows one example that erroneously detects the cytoplasm due to the incorrect a* channel enhancement result. Actually, the performance of a majority of our H&E stained images resembles the examples in Figs. 6 and 8. The proposed nucleus segmentation method was designed with real-time application of cervical cancer screening in mind. Hence the LAGC and ellipse fitting only process a part of nuclei (as illustrated in Figs. 4 and 5, respectively), and thus it may reduce the sensitivity to segmentation of abnormal nuclei. Moreover, current nucleus segmentation still may not be able to accurately detect the boundary of some nuclei with poor staining and contrast (Fig. 11(a)), and may erroneously detect boundary of some nuclei in atrophic cells (Fig. 11(b)).




11

Fig. 10. One example of erroneous segmentation of some cytoplasm due to intensively bright illumination. (a) Original image, (b) ground truth, (c) a* channel enhancement, (d) cytoplasm segmentation result. Arrow points to the bright red region and therefore erroneous segmented cytoplasm.

Fig. 11. Two examples of nuclei with poor staining and contrast (a) and in atrophic cell (b). In each graph, the first through the third columns illustrate the original image, nucleus segmentation by LAGC, and ground truth.

5. Conclusion

Acknowledgements

To the best of our knowledge, this paper is the first attempt to achieve automated segmentation of both healthy and abnormal cervical cells in a FOV, given existing techniques more often designed for healthy cells and ideal imaging conditions. The proposed global GC method allows for cytoplasm delineation even when image histograms present non-bimodal distribution. The proposed LAGC approach and concave points-based method enable the nucleus segmentation in images with pathology and overlapping cells. The method was tested on 21 cervical cell images with realistic, non-ideal imaging condition and pathology (abnormal nuclei). The experimental results yielded a cytoplasm segmentation accuracy of DSC = 93% and an abnormal nuclei binarization accuracy of F-measure = 88.4%. The proposed method can complete the cytoplasm and nuclei segmentation tasks in a FOV in 1.9 s, making the proposed method highly practical in most clinical settings.

The authors would like to thank Jingli Li, Shuangming Zheng, and Sheng Tang at the Shenzhen Microprofit Electronic Co. Ltd., for their participation in this work; Milan Sonka, Andreas Wahle, Li Zhang, and Junjie Bai from the University of Iowa, for their contributions to the manuscript revision and English writing. This work was supported by the Key Program of National Natural Science Foundation of China (Grant Number: 61031003), and the National Natural Science Foundation of China (Grant Numbers: 61372006, 61271108). References [1] Saslow D, Solomon D, Lawson HW, Killackey M, Kulasingam SL, Cain J, Garcia FA, Moriarty AT, Waxman AG, Wilbur DC, Wentzensen N, Downs Jr LS, Spitzer M, Moscicki AB, Franco EL, Stoler MH, Schiffman M, Castle PE, Myers ER. The ACS-ASCCP-ASCP Cervical Cancer Guideline Committee. American cancer society, American Society for Colposcopy and Cervical Pathology, and American



[2] [3]

[4] [5] [6]

[7]

[8] [9]

[10]

[11] [12] [13] [14] [15]

[16]

[17]

[18] [19] [20]

[21]

[22] [23]


Society for Clinical Pathology Screening Guidelines for the Prevention and early detection of cervical cancer. CA Cancer J Clin 2012;62:147–72. Birdsong GG. Automated screening of cervical cytology specimens. Hum Pathol 1996;27:468–81. Kitchener HC, Blanks R, Dunn G, Gunn L, Desai M, Albrow R, Mather J, Rana DN, Cubie H, Moore C, Legood R, Gray A, Moss S. Automation-assisted versus manual reading of cervical cytology (MAVARIC): a randomised controlled trial. Lancet Oncol 2011;12:56–60. Tsai MH, Chan YK, Lin ZZ, Yang-Mao SF, Huang PC. Nucleus and cytoplast contour detector of cervical smear image. Pattern Recogn Lett 2008;29:1441–53. Li K, Lu Z, Liu W, Yin J. Cytoplasm and nucleus segmentation in cervical smear images using Radiating GVF Snake. Pattern Recogn 2012;45:1255–64. Yang-Mao SF, Chan YK, Chu YP. Edge enhancement nucleus and cytoplast contour detector of cervical smear images. IEEE Trans Syst Man Cy B 2008;38:353–66. Bergmeir C, García Silvente M, Benítez JM. Segmentation of cervical cell nuclei in high-resolution microscopic images: a new algorithm and a web-based software framework. Comput Method Prog Biol 2012;107:497–512. Gençtav A, Aksoy S, Önder S. Unsupervised segmentation and classification of cervical cell images. Pattern Recogn 2012;45:4151–68. Harandi NM, Sadri S, Moghaddam NA, Amirfattahi R. An automated method for segmentation of epithelial cervical cells in images of ThinPrep. J Med Syst 2010;34:1043–58. Lu Z, Carneiro G, Bradley AP. Automated nucleus and cytoplasm segmentation of overlapping cervical cells. Med Image Comput Comput Assist Interv 2013:452–60. Jantzen J, Norup J, Dounias G, Bjerregaard B. Pap-smear benchmark data for pattern classification. Proc Nat Insp Smart Inf Syst 2005. Hu M, Ping X, Ding Y. Automated cell nucleus segmentation using improved snake. Int Conf Image Proc 2004:2737–40. Bamford P, Lovell B. Unsupervised cell nucleus segmentation with active contours. Signal Process 1998;71:203–13. Wu HS, Barba J, Gil J. A parametric fitting algorithm for segmentation of cell images. IEEE Trans Bio-Med Eng 1998;45:400–7. Lee SJ, Wilhelm PS, Meyer MG, Bannister WR, Kuan CL, Ortyn WE, Nelson LA, Frost KL, Hayenga JW. Cytological slide scoring apparatus. Neo Path, Inc., U.S. Pub. No. 5933519, August 3, 1999. Plissiti ME, Nikou C, Charchanti A. Automated detection of cell nuclei in Pap smear images using morphological reconstruction and clustering. IEEE Trans Inf Technol B 2011;15:233–41. Plissiti ME, Nikou C, Charchanti A. Combining shape, texture and intensity features for cell nuclei extraction in Pap smear images. Pattern Recogn Lett 2011;32:838–53. Jung C, Kim C, Chae SW, Oh S. Unsupervised segmentation of overlapped nuclei using Bayesian classification. IEEE Trans Bio-Med Eng 2010;5:2825–32. Plissiti ME, Nikou C. Overlapping cell nuclei segmentation using a spatially adaptive active physical model. IEEE Trans Image Process 2012;21:4568–80. Chang CW, Lin MY, Harn HJ, Harn YC, Chen CH, Tsai KH, Hwang CH. Automatic segmentation of abnormal cell nuclei from microscopic image analysis for cervical cancer screening. IEEE Int Conf Nano/Molecular Med Eng 2009:18–21. Nielsen B, Albregtsen F, Danielsen HE. Automatic segmentation of cell nuclei in Feulgen-stained histological sections of prostate cancer and quantitative evaluation of segmentation results. Cytom Part A 2012;81:588–601. Boykov Y, Veksler O, Zabih R. Fast approximate energy minimization via graph cuts. IEEE Trans Pattern Anal Mach Intel 2001;23:1222–39. Al-Kofahi Y, Lassoued W, Lee W, Roysam B. Improved automatic detection and segmentation of cell nuclei in histopathology images. IEEE Trans Bio-Med Eng 2010;57:841–52.

[24] Lou X, Koehe U, Wittbrodt J, Hamprecht FA. Learning to segment dense cell nuclei with shape prior. IEEE Conf Comput Vis Pattern Recogn 2012: 1012–8. [25] Chang H, Han J, Borowsky A, Loss L, Gray JW, Spellman PT, Parvin B. Invariant delineation of nuclear architecture in glioblastoma multiforme for clinical and molecular association. IEEE Trans Med Image 2013;32: 670–82. [26] Gonzalez RC, Woods RE. Digital image processing. New Jersey: Pearson Education, Inc; 2008. p. 954. [27] Kong H, Gurcan M, Belkacem-Boussaid K. Partitioning histopathological images: an integrated framework for supervised color-texture segmentation and cell splitting. IEEE Trans Med Image 2011;30:1661–77. [28] Zhang L, Chen SP, Wang TF, Chen Y, Liu SX, Li MH. A practical segmentation method for automated screening of cervical cytology. Proc IEEE Int Conf Intell Comput Bio-Med Inst 2011:140–3. [29] Otsu NA. Threshold selection method from gray-level histograms. IEEE Trans Syst Man Cy 1979;9:62–6. [30] Boykov Y, Kolmogorov V. An experimental comparison of min-cut/max-flow algorithms for energy minimization in vision. IEEE Trans Pattern Anal Mach Intel 2004;26:1124–37. [31] Bishop C. Pattern recognition and machine learning. New York: Springer; 2007. p. 738. [32] Sahoo PK, Soltani S, Wong AKC, Chen Y. A survey of thresholding techniques. Comput Vis Graph Image Process 1988;41:233–60. [33] Sauvola J, Pietikainen M. Adaptive document binarization. Pattern Recogn 2000;33:225–36. [34] Shafait F, Keysers D, Breuel TM. Efficient implementation of local adaptive thresholding techniques using integral images. Proc SPIE, Doc Recogn Retr 2008;XV:681510–5. [35] Van der Laak JA, Siebers AG, Cuijpers VM, Pahlplatz MM, de Wilde PC, Hanselaar AG. Automated identification of diploid reference cells in cervical smears using image analysis. Cytom 2002;47:256–64. [36] Trier BD, Jain AK. Goal-directed evaluation of binarization methods. IEEE Trans Pattern Anal Mach Intel 1995;17:1191–201. [37] Viola P, Jones MJ. Robust real-time face detection. Int J Comput Vision 2004;57:137–54. [38] Haralick RM. A measure of circularity of digital figures. IEEE Trans Syst Man Cybern 1974:394–6. SMC-4. [39] Loy G, Zelinsky A. Fast radial symmetry for detecting points of interest. IEEE Trans Pattern Anal Mach Intel 2003;25:959–73. [40] Maksem JA, Finnemore M, Belsheim BL, Roose EB, Makkapati SR, Eatwell L, Weidmann J. Manual method for liquid-based cytology: a demonstration using 1,000 gynecological cytologies collected directly to vial and prepared by a smear-slide technique. Diagn Cytopathol 2001;25:334–8. [41] Dice LR. Measures of the amount of ecologic association between species. Ecology 1945;26:297–302. [42] Comaniciu D, Meer P. Mean-shift: a robust approach toward feature space analysis. IEEE Trans Pattern Anal Mach Intel 2002;24:603–19. [43] Chan TF, Vese LA. Active contours without edges. IEEE Trans Image Process 2001;10:266–77. [44] Cox DR, Hinkley DV. Theoretical statistics. Cambridge, MA: Oxford University; 1974. p. 528. [45] Bland JM, Altman DG. Statistical methods for assessing agreement between two methods of clinical measurement. Lancet 1986;1:307–10. [46] Marinakis Y, Dounias G, Jantzen J. Pap smear diagnosis using a hybrid intelligent scheme focusing on genetic algorithm based feature selection and nearest neighbour classification. Comput Biol Med 2009;39:69–78.


Accurate Segmentation of Cervical Cytoplasm and Nuclei Based on Multiscale Convolutional Network and Graph Partitioning.

SEGMENTATION OF MYOCARDIUM USING DEFORMABLE REGIONS AND GRAPH CUTS.

Prostate MRI segmentation using learned semantic knowledge and graph cuts.

Cytoplasm segmentation on cervical cell images using graph cut-based approach.

A multiatlas segmentation using graph cuts with applications to liver segmentation in CT scans.

A deep learning based framework for accurate segmentation of cervical cytoplasm and nuclei.

Abnormal cervical cytology.

Abnormal cervical cytology.

Bayesian segmentation of atrium wall using globally-optimal graph cuts on 3D meshes.

Automatic Liver Segmentation on Volumetric CT Images Using Supervoxel-Based Graph Cuts.

The use of atlas registration and graph cuts for prostate segmentation in magnetic resonance images.

Segmentation of anterior cruciate ligament in knee MR images using graph cuts with patient-specific shape constraints and label refinement.

Segmentation of choroidal boundary in enhanced depth imaging OCTs using a multiresolution texture based modeling in graph cuts.

Semi-automatic breast ultrasound image segmentation based on mean shift and graph cuts.

Susceptibility artefact correction using dynamic graph cuts: application to neurosurgery.

Segmentation of MR image using local and global region based geodesic model.

Drosophila Eye Nuclei Segmentation Based on Graph Cut and Convex Shape Prior.

Colposcopic evaluation of patients with abnormal cervical cytology.

Immunoglobulin concentrations in cervical mucus in patients with normal and abnormal cervical cytology.

Colposcopic management of abnormal cervical cytology during pregnancy.

Abnormal cervical cytology leading to the diagnosis of gastric cancer.

Evaluation of abnormal cervical cytology during pregnancy with colposcopy.

Untypable human papillomavirus infection and risk of cervical intraepithelial neoplasia among women with abnormal cervical cytology.

Abnormal cervical cytology in the teen-ager: a continuing problem.