Semantic Segmentation of Colon Glands with Deep Convolutional Neural Networks and Total Variation Segmentation

Introduction

The variability of glandular structures in biological tissue poses a challenge to automated analysis of histopathology slides. It has become a key requirement to quantitative morphology assessment and supporting cancer grading. Considering non-pathological cases only, automated segmentation algorithms must already be able to deal with significant variability in shape, size, location, texture and staining of glands. Moreover, in pathological cases gland objects can tremendously differ from non-pathological and benign glands, which further exacerbates finding a general solution to the segmentation problem.

Previous work on gland segmentation in colon tissue has used graphical models or textural features . Others worked on segmentation in prostatic cancer tissue using an integrated low-, high-level and contextual segmentation model , probabilistic Markov models , k-means clustering and region growing , spatial association of nuclei to gland lumen . The reader is referred to the work of Sirinukunwattana et al. for a more detailed description of work related to glandular structure segmentation. Deep learning methods, especially convolutional neural networks (CNNs) , have found applications in biomedical image analysis for different tasks: semantic segmentation , mitosis detection and classification , and blood cell counting .

In this work, we propose a learning-based strategy to semantically segment glands in the Warwick-QU dataset, presented at the GlaS@MICCAI2015 challengehttp://www2.warwick.ac.uk/fac/sci/dcs/research/combi/research/bic/glascontest/. It contains 161 annotated images of benign and malignant colorectal adenocarcinoma, stained with Hematoxylin-Eosin (H&E) and scanned at $20\times$ magnification. Fig. 1 shows some example images and their ground truth annotation. In each image, individual objects are annotated with the same label, illustrated by the different colors. To the challenge participants, information on whether an image shows benign or malignant tissue is only available in the training dataset. Three datasets were released during the contest and the total number of non-overlapping images (benign/malignant) in the training set, test set A and test set B is $85(37/48)$, $60(33/27)$, and $16(12/4)$, respectively. These datasets further contained $795$, $666$, and $91$ individual glands.

The contributions of our work are twofold: (i) we present a novel deep learning scheme to generate classifier predictions for malignant and benign object and background pixels accompanied by a dedicated gland-separating refinement classifier that is able to distinguish touching objects, which pose a challenge for later segmentation. (ii) We use these classification results as the input for a simple, yet effective, globally optimal figure-ground segmentation approach based on a convex geodesic active contour formulation that regularizes the classifier predictions according to a minimal contour-length principle. Both technological contributions are described in section 2, while the subsequent sections show and discuss the results of our novel approach applied to the datasets of the GlaS@MICCAI2015 challenge.

Methods

We present a segmentation method for Hematoxylin-Eosin (H&E) stained histopathological sections that proceeds in three steps: The raw RGB images are preprocessed to extract a robust representation of the tissue structure. Subsequently, two classifiers are trained to predict glands (Object-Net) and gland-separating structures (Separator-Net) from the image. Finally, the outputs of the classifiers are combined and a figure-ground segmentation based on weighted total variation is used to produce the segmentation result.

Prior to classification, the RGB images are preprocessed as shown in Fig. 2. A standard color deconvolution is performed for the specific H&E staining used in the provided datasetWe used the H&E 2 setting in the implementation of G. Landini, available in Fiji .. It separates tissue components according to their staining, emphasizes the structure and inherently performs data whitening. The first (red) channel of the deconvolved RGB image contains most of the tissue structure information, so the other channels can be omitted. In order to account for different staining contrasts and lighting conditions during image acquisition, contrast limited adaptive histogram equalization (CLAHE) is applied.

2 Learning Pixel Classifiers

Given the large variability of both benign and malignant tissue in the Warwick-QU dataset, we opted for CNNs due to their recently shown convincing performance in pixelwise classification of histopathology images and to learn a rich set of features directly from images.

The general architecture of both CNNs is motivated by a classical LeNet-5 architecture and consists of $K=7$ ($k=1,\ldots,K$) layers: four convolutional layers (Conv$k$) for feature learning and three fully connected (FC$k$) layers as feature classifier, see Fig. 3. The rectified linear unit (ReLU) nonlinearity ($f(x)=\max(0,x)$) is used as the activation function throughout all layers in the networks. All convolutional layers consist of a set of learnable square 2D filters with pixel stride 1, followed by ReLU activation. Subsampling (max-pooling) layers (Sub$k$, $2\times 2$), accounting for translation invariance, are used after the first three convolutional layers and are counted as part of the convolutional layer. The final pixelwise classification of an input image is obtained by sliding a window over the image, and classifying the center pixel of that window.

For training minibatch stochastic gradient descent (MBSGD) with momentum, weight decay, and dropout regularization is used to minimize a negative log-likelihood loss function.

The goal of the Object-Net is to predict the probability of a pixel belonging to a gland or background. One could now define a binary classification problem, but malignant and benign tissue express express unique features, which are not found in the other tissue type, and which can thus complicate the learning problem. We therefore formulate an alternative, four-class classification problem, in which we distinguish ($L=4$, with $l=0,...,L-1$): background benign ($C_{0}$), gland benign ($C_{1}$), background malignant ($C_{2}$), and gland malignant ($C_{3}$). In order to do that it is necessary to transform the provided ground truth labels to reflect benignity and malignancy as well. The annotation images are binarized and a new label is assigned to pixels belonging to each class $C_{l}$, see Fig. 4.

The input to the CNN is an image patch $I(\mathbf{x})$ of size $101\times 101$ pixels, centered at an image location $\mathbf{x}=(u,v)^{\top}$, where $\mathbf{x}\in\Omega$ and $\Omega$ denotes the image domain. A given patch $I(\mathbf{x})$ is convolved with 80 filters ($11\times 11$) in the first convolutional layer, in the second layer with 96 filters ($7\times 7$), in the third layer with 128 filters ($5\times 5$), and in the last layer with 160 filters ($3\times 3$), see Fig. 3(a). The three subsequent fully connected layers FC5-FC7 of the classifier contain 1024, 512, and four output units, respectively. The output of FC7 is fed into a softmax function, producing the center pixel’s probability distribution over the labels. The probability for each class $l$ is stored in a corresponding map $I_{C_{l}}(\mathbf{x})$.

2.2 Separator-Net: Classifying Gland-separating Structures

Initial experiments have shown that taking pixelwise predictions only from the Object-Net were insufficient in order to separate very close gland objects. Hence, a second CNN, the Separator-Net, is trained to predict structures in the image that are separating such objects. This learning problem is formulated as binary classification task.

As depicted in Fig. 3(b), the CNN structure is similar to the Object-Net: a given input image patch $I(\mathbf{x})$ of size $101\times 101$ pixels is convolved with 64 filters ($9\times 9$) in the first convolutional layer, in the second layer with 96 filters ($7\times 7$), in the third layer with 128 filters ($5\times 5$), and in the last layer with 160 filters ($3\times 3$). The three subsequent fully connected layers FC5-FC7 of the classifier contain 1024, 512, and two output units, respectively. The output of the last layer (FC7) is fed into a softmax function to produce the probability distribution over the labels for the center pixel. The probability for a pixel $\mathbf{x}$ belonging to a gland-separating structure is stored in the corresponding probability map $S(\mathbf{x})$.

2.3 Refining CNN Outputs

Once all probability maps have been obtained, the Object-Net predictions $I_{C_{l}}(\mathbf{x})$ are refined with the Separator-Net predictions $S(\mathbf{x})$ to emphasize the gland borders and prevent merging of close objects. The subsequent figure-ground segmentation algorithm requires a single foreground and background map to produce the final segmentation result, so outputs are combined as follows.

The foreground probability map $p_{fg}$ is constructed by

where $\rho\in\left[0,1\right]$ controls the influence of the refinements done by the separator predictions. Similarly, evaluating Eq. (2) produces the background probability map:

3 Total Variation Segmentation

To generate a final segmentation, the following continuous non-smooth energy functional $E_{seg}(u)$ is minimized:

where $\nabla I(\mathbf{x})$ is the gradient of the input image, thus attracting the segmentation towards large gradients. The second term in Eq. (3) is the data term with $w$ describing a weighting map. The values in $w$ have to be chosen negative if $u$ should be foreground and positive if $u$ should be background. If values in $w$ are set to zero, the pure weighted TV energy is minimized seeking for a minimal contour length segmentation. We use the refined outputs from the previous classification step (Eqs. (1) and (2)) and introduce a threshold $\tau$ to ensure a minimum class confidence in a map $p$:

The weighting map $w$ is derived by applying the logit transformation:

The regularization parameter $\lambda$ defines the trade-off between our data term and the weighted TV semi-norm. The stated convex problem in Eq. (3) can be solved for its global optimum efficiently using the primal-dual algorithm , which can be implemented very efficiently using NVidia CUDA, thus making use of the parallel computing power of recent GPUs. As the segmentation $u$ is continuous, the final segmentation is achieved by thresholding $u$ with a value of $0.5$. We optimize the free parameters $\alpha$, $\beta$ and $\lambda$ by performing a grid search in a suitable range of these values ($\alpha\in[0.5,15]$, $\beta\in[0.35,0.95]$ and $\lambda\in[0.01,10]$), where all 85 annotated training images are used to tune these parameters based on the Dice coefficient.

4 Implementation Details

For the sake of execution speed when using a sliding window approach, the images are rescaled to half resolution prior to classification and upsampled with bilinear interpolation to their original size afterwards. The size of the input patch $I(\mathbf{x})$ is chosen to be $101\times 101$ pixels, such that sufficient contextual information is available to classify the center pixel.

The majority of training images (79) have a size of $775\times 522$ pixels, and resizing reduces them to $387\times 261$ pixels. If we just considered the valid part without border extension for sampling the patches for the training dataset, we would actually lose approximately $46\%$ of the labeled pixels when using a patch size of $101\times 101$ pixels. On the other hand, we would introduce a significant number of boundary artifacts by artificially extending the border to make use of all labeled pixels. Fortunately, most images are tiles of a bigger image and can thus be stitched seamlessly to obtain a total of 19 imagesIn one case, stitching was not possible, since only 3 tiles were available. These 3 tiles, and the remaining 6 images, that were not part of a bigger image, were treated as individual images. (Fig. 5), where we can sample enough patches without heavily relying on artificial border extension.

In principle, we pursued the same sampling strategy for the Separator-Net, but were required to create the ground truth labels manually. We annotated all pixels that belong to a structure very close to two or more gland borders. The green lines in Fig. 5 illustrate the additional manual annotation of the separating structures. Due to the low number of foreground samples when compared to the Object-Net, the number of foreground samples for the Separator-Net was artificially increased by exploiting the problem’s requirement for rotation-invariance and adding nine additional rotated versions of the patch, i.e. every $36^{\circ}$.

4.2 CNN Training

Both CNNs were trained on a balanced training set of $125,000$ image patches per class. Patches in the training sets were sampled at random from the available pool of training images. Training and test sets reflect approximately the same distribution of samples over images. The size of the minibatches in the MBSGD was set to $200$ samples and the networks were trained until the stopping criterion was met: no further improvement of the error rate on a held-out test set over 20 epochs. We set the initial learning rate $\eta_{0}=0.0025$, with a linear decay saturating at $0.2\eta_{0}$ after 100 epochs. For all layers, a weight decay was chosen to be $0.005$ and the dropout rate was set to $0.5$. We used an adaptive momentum term starting at $0.8$ and increasing to $0.99$ after 50 epochs, such that with progressing training the updates are influenced by a larger number of samples than at the beginning.

Fig. 6 shows the classification error rate as a function of the training duration in epochs. Each class was represented with $5,000$ samples in the test set for the Object-Net, and $10,000$ for the Separator-Net, respectively. The training error is actually estimated on a fixed subset of the training data ($20,000$ samples), to get an intuition when overfitting starts. The Object-Net achieves the best performance after $43$ epochs, with a minimum training error of $0.0475$ and a minimum test error of $0.0492$. Training of the Separator-Net continued until the lowest training error of $0.0231$ and test error of $0.0624$ was reached after $119$ epochs. Fig. 7 shows the learned filters of the first convolutional layer in both networks. The CNN models were implemented in Pylearn2 , a machine learning library built on top of Theano .

Results

The grid search resulted in $\alpha=10$, $\beta=0.95$ and $\lambda=0.1$ as parameters optimizing the TV segmentation based on the Dice score. The confidence threshold for foreground and background was determined empirically and fixed to $\tau=0.65$. Separator predictions were fully considered for refining the Object-Net predictions ($\rho=1$).

In Table 1, we report performance metricsThe evaluation scripts were kindly provided by the contest organizers and are available from http://www2.warwick.ac.uk/fac/sci/dcs/research/combi/research/bic/glascontest/ evaluation/. for detection (precision, recall, F1-score), segmentation (object-level Dice), and shape (Hausdorff distance) on the training set, as well as test set A and B as mean and standard deviation (SD). Blobs with an area less than $500$ pixels were removed and all remaining blobs were labeled with unique identifiers before computing the measures.

Compared to using predictions only from the Object-Net, the segmentation performance improved with separator refinement. Malignant cases are harder to segment due to their irregular shape and pathological variations in the tissue. Fig. 8 illustrates some qualitative example segmentation results on the training data set, Fig. 9 and Fig. 10 show results on test set A and B, respectively.

The average total runtime for segmenting a $577\times 522$ image is 5 minutes using an NVidia GeForce Titan Black 6GB GPU.

2 Benignity and Malignancy Classification

In the proposed approach, the Object-Net inherently learns a discrimination of benign ($c=0$) and malignant ($c=1$) tissue, since the labels for benign and malignant are available in the training dataset and we defined a four-class classification problem. Instead of combining the probability maps for glands and background as done for segmentation, we combine the maps for benignity and malignancy. Subsequently, the average probabilities for a benign case can be computed as

where $|\Omega|$ is the number of pixels in the image domain $\Omega$. The maximum of both values finally indicates the prediction:

We evaluated the classification performance for benign and malignant tissue on the two test sets A and B and achieved an accuracy of $98.33\%$ and $93.75\%$. The average (SD) decision confidence in test set A was $0.84(0.13)$ for benign and $0.81(0.11)$ for malignant, and in test set B $0.74(0.11)$ and $0.86(0.15)$, respectively.

Discussion and Conclusions

This paper presented a method to segment glands in H&E stained histopathological images of colorectal cancer using deep convolutional neural networks and total variation segmentation. As our main contribution, we showed that segmentation results can be greatly improved when the predictions of the Object-Net are refined with the learned gland-separating structures of the Separator-Net. Adding the separators does not only regulate the trade-off between precision and recall, but generally improves the performance scores for detection (F1-score), segmentation (Dice) and shape (Hausdorff). The final ranking as well as the test set performance results of other algorithms participating in this challenge are available online at the contest websitehttp://www2.warwick.ac.uk/fac/sci/dcs/research/combi/research/bic/glascontest/ results/, which is continuously being updated by algorithms from new participating groups.

Our approach inherently allows to very accurately discriminate benign and malignant cases, because the Object-Net was trained on labels for both cases. The average confidence for a decision towards benignity and malignancy is acceptable. Nevertheless, we cannot distinguish more detailed histologic grades among these cases, since there was no information (e.g. high- or low-grade) available in addition to the segmentation ground truth.

Acknowledgements

The authors are grateful to the organizers of the GlaS@MICCAI2015 challenge for providing (i) the Warwick-QU image dataset, and (ii) the MATLAB evaluation scripts for computing performance measures that are comparable among the participating teams. Further thanks goes to Julien Martel for fruitful discussions in early phases of this challenge.

References