Abdominal multi-organ segmentation with organ-attention networks and statistical fusion

Abstract

Accurate and robust segmentation of abdominal organs on CT is essential for many clinical applications such as computer-aided diagnosis and computer-aided surgery. But this task is challenging due to the weak boundaries of organs, the complexity of the background, and the variable sizes of different organs. To address these challenges, we introduce a novel framework for multi-organ segmentation of abdominal regions by using organ-attention networks with reverse connections (OAN-RCs) which are applied to 2D views, of the 3D CT volume, and output estimates which are combined by statistical fusion exploiting structural similarity. More specifically, OAN is a two-stage deep convolutional network, where deep network features from the first stage are combined with the original image, in a second stage, to reduce the complex background and enhance the discriminative information for the target organs. Intuitively, OAN reduces the effect of the complex background by focusing attention so that each organ only needs to be discriminated from its local background. RCs are added to the first stage to give the lower layers more semantic information thereby enabling them to adapt to the sizes of different organs. Our networks are trained on 2D views (slices) enabling us to use holistic information and allowing efficient computation (compared to using 3D patches). To compensate for the limited cross-sectional information of the original 3D volumetric CT, e.g., the connectivity between neighbor slices, multi-sectional images are reconstructed from the three different 2D view directions. Then we combine the segmentation results from the different views using statistical fusion, with a novel term relating the structural similarity of the 2D views to the original 3D structure. To train the network and evaluate results, $13$ structures were manually annotated by four human raters and confirmed by a senior expert on $236$ normal cases. We tested our algorithm by 4-fold cross-validation and computed Dice-Sørensen similarity coefficients (DSC) and surface distances for evaluating our estimates of the $13$ structures. Our experiments show that the proposed approach gives strong results and outperforms 2D- and 3D-patch based state-of-the-art methods in terms of DSC and mean surface distances.

Introduction

Segmentation of the internal structures, like body organs, in medical images is an essential task for many clinical applications such as computer-aided diagnosis (CAD), computer-aided surgery (CAS) and radiation therapy (RT). However, despite intensive studies of automatic or semi-automatic segmentation methods, there remain challenges which need to be overcome before these methods can be applied to clinical environments. In particular, detailed abdominal organ segmentation on CT is a challenging task both for manual human annotation and for automatic segmentation algorithms for various reasons including the morphological complexity of the structures, the large variations between inter- and intra-subjects, and image characteristics such as low contrast of soft tissues.

Early studies of abdominal organ segmentation focused on specific single organs, for example relatively large isolated structures such as the liver or critical structures such as blood vessels . However, most of the algorithms were based on specific features of the target organ, and so extensibility to the simultaneous segmentation of multiple organs was limited. For multi-organ segmentation, atlas-based approaches were adopted for many applications . The general framework of atlas-based segmentations is to deformably register selected atlas images with segmented structures to the target image. Critical issues for this approach, which affect performance accuracy, include proper atlas selection, accurate deformable image registration, and label fusion. In particular, for the abdominal region, inter-subject variations are relatively large compared with other parts of the body (e.g., the brain) so the segmentation results are dependent on deformable registration between inter-subjects from the limited set of atlases, which is a challenging problem that critically affects the final accuracies. In addition, computational time is strongly dependent on the number of atlases. Therefore, selection of the proper number and types of atlases is a critical factor for both of the accuracy and efficiency.

Recently, learning-based approaches exploiting large datasets have been applied to the segmentation of medical images . In particular, deep convolutional neural networks (CNN) have been very successful . Targets include regions in the brain , chest , and abdomen . The performance results of CNNs for organs (and even tumors) reach, or outperform, alternative state-of-the-art methods. Unlike multi-atlas-based approaches, deep networks do not require selecting a specific atlas or require deformable registration from training sets to a target image. In this study, we apply deep network approaches to abdominal organ segmentation.

Most studies based on deep networks, however, focused on a single structure segmentation, particularly for abdominal regions, and there are few studies of multi-organ segmentation partly due to technical challenges discussed later. We note that fully convolutional networks (FCNs) have been generally accepted for organ segmentations on CT scans partly because they give state-of-the-art performance for semantic segmentation of natural images . But there are three major characteristics of abdominal CT which we must address in order to obtain strong performance on multi-organ segmentation.

Firstly, many abdominal organs have weak boundaries between spatially adjacent structures on CT, e.g. between the head of the pancreas and the duodenum. In addition, the entire CT volume includes a large variety of different complex structures. Morphological and topological complexity includes anatomically connected structures such as the gastrointestinal (GI) track (stomach, duodenum, small bowel and colon) and vascular structures. The correct anatomical borders between connected structures may not be always visible in CT, especially in sectional images (i.e., 2D slices), and may be indicated only by subtle texture and shape change, which causes uncertainty even for human experts. This makes it hard for deep networks to distinguish the target organs from the complex background.

Secondly, there are large variations in the relative sizes of different target organs, e.g. the liver compared to the gallbladder. This causes problems when applying deep networks to multi-organ segmentation because lower layers typically lack semantic information when segmenting small structures. The same problem has been observed in semantic segmentation of natural images where the segmentation performance on small regions is typically much worse than on large regions, motivating the need to introduce mechanisms which attend to the scale .

Thirdly, although CT scans are high-resolution three-dimensional volumes, most current deep network methods were designed for 2D images. To overcome the limitations of using 2D CNNs for 3D images, Setio et al. used multiple 2D patches reconstructed from $9$ different directions around the target region for the task of pulmonary nodule detection. Zhuang et al. used 2D axial, coronal, and sagittal slices for pancreas detection at the coarse level and also for segmentation at the finer level. More recently, there are studies which use 3D deep networks . These, however, are not networks that act on the entire 3D CT volume but instead are local patch-based approaches (due to complex challenges of 3D deep networks discussed later in this paragraph). To address the problems caused by restricting to image patches, used a hierarchical approach with multi-resolutions, which reduces the dimension of the whole volume for initial detection and focuses on smaller regions at the finer resolution. But this strategy is best suited to a single target structure. Roth et al. applied a bigger patch size to deal with the whole dense pancreatic volume, but this was also for single pancreas segmentation and hard to extend to the whole abdominal region. In general, 3D deep networks face far greater complex challenges than 2D deep networks. Both approaches rely heavily on graphics processing units (GPUs) but these GPUs have limited memory size which makes it difficult when dealing with full 3D CT volumes compared to 2D CT slices (which require much less memory). In addition, 3D deep networks typically require many more parameters than 2D deep networks and hence require much more training data, unless they are restricted to patches. But there is limited training data for abdominal CT images, because annotating them is challenging and requires expert human radiologists, which makes it particularly difficult to apply 3D deep networks to abdominal multi-organ segmentation. We have, however, implemented a 3D patch based approach for comparison.

To deal with the technical difficulties for abdominal multi-organ segmentation on CT, we introduce a novel framework of an organ-attention 2D deep networks with reverse connections (OAN-RC) followed by statistical fusion to combine the information from the three different views exploiting structural similarity using local isotropic 3D patches. OAN is a two-stage deep network, which computes an organ-attention map (OAM) from typical probability map of labels for input images in the first stage and combines OAM to the original input image for the second stage. This two-stage strategy effectively reduces the complexity of the background while enhancing the discriminative information of target structures (by concentrating attention close to the target structures). By training OAM with additional deep network, uncertainties and errors from the first stage are adjusted and the fidelity of the final probability map is improved. In this procedure, we apply reverse connections to the first stage so that we can localize organ information at different scales by assisting the lower layers with semantic information.

More specifically, we apply OAN-RC to each sectional slice, which is an extreme form of anisotropic local patches but include the whole semantic (i.e. volume) information from one viewing direction. This yields segmentation information from separate sets of multi-sectional images (axial, coronal, and sagittal planes in this study similarly to most of medical image platforms for 2D visualization). We statistically fuse the three sources of information using local isotropic 3D patches based on direction-dependent local structural similarity. The basic fusion framework uses expectation-maximization (EM) similar to . But, unlike typical statistical fusion methods used for atlas-based segmentation, the input volumes and the target volumes for segmentation in our problem are the same. But different structures and texture patterns, from different viewing directions, will often generate nonidentical segmentations in 3D. Our strategy is to exploit structural similarity by computing a direction-dependent local property at each voxel. This models the structural similarity from the 2D images to the original 3D structure (in the 3D volume) by local weights. This structural statistical fusion improves our overall performance by combining the information from the three different views in a principled manner and also imposing local structure.

Figure 1 describes the graphical concept of our framework. Our proposed algorithm was tested on $236$ abdominal CT scans of normal cases collected as a part of FELIX project for pancreatic cancer research . By experiments, our method showed robust and high fidelities to the ground-truth for all target structures with smooth boundaries. It outperformed 3D patch-based algorithms as well as 2D-based in terms of DICE-similarity coefficient and average surface distance with memory and computational efficiency.

Organ-Attention Networks with Reverse Connections

We first introduce the OAN, which is composed of two jointly optimized stages. The first stage (stage-I) transforms the organ segmentation probability map to provide spatial attention to the second stage (stage-II), so that the segmentation network trained in stage-II is more discriminative for segmenting organs (because it only has to deal with local context). To assist the lower layers in stage-I with more semantic information, we employ reverse connections (Sec. 2.2), which pass semantic information down from high layers to low layers. The OAN is trained in an end-to-end fashion to enhance the learning ability of all stages.

The input images to our OAN are reconstructed 2D slices from axial, sagittal and coronal directions. Based on the normal vector directions of the sagittal ($X$), coronal ($Y$) and axial ($Z$) planes, we denote the 2D images by $\mathbf{I}_{i}^{X}$, $\mathbf{I}_{j}^{Y}$ and $\mathbf{I}_{k}^{Z}$ respectively, where $i=1,\ldots,n_{x},~{}j=1,\ldots,n_{y},~{}k=1,\ldots,n_{z}$ and $n_{x},~{}n_{y},~{}n_{z}$ are the numbers of slices for the three directions, respectively, and $\bigcup_{i}\mathbf{I}^{X}_{i}=\bigcup_{j}\mathbf{I}^{Y}_{j}=\bigcup_{k}\mathbf{I}^{Z}_{k}=V$. Following the work of , we train an individual OAN for each direction.

where $\mathbf{1}(\cdot)$ is an indicator function.

Using a preliminary organ segmentation map to guide the computation of a better organ segmentation can be thought as employing an attentional mechanism. Towards this end, we propose an organ-attention module by

where $*$ denotes the convolution operator, $\mathbf{W}$ indicates the convolutional filters, and $\mathbf{b}$ is the bias. (2) embeds cross-organ information into a single organ-attention map, $\mathbf{Q}$, which learns discriminative spatial attention for different organs automatically. By combining $\mathbf{Q}$ with the original input $\mathbf{I}$, we get an image which emphasizes each organ by

where $\star$ is the element-wise product operator. We apply $\mathbf{I}^{(2)}$ to the input of stage-II, and the probability of stage-II then becomes $\mathbf{P}^{(2)}=f(\mathbf{I}^{(2)};\bm{\Theta}^{(2)})$.

In order to drive stage-II to focus on organ regions without needing to deal with complicated non-local background, we define a selection function, $\mathbf{1}(\mathbf{P}^{(1)}_{0}\leqslant\rho)$ where $\mathbf{P}^{(1)}_{0}=\{{p}_{i,0}^{(1)}\}_{i=1,...,H\times W}$ is the probability map provided by stage-I. In stage-II, we only accept the region if $p_{i,0}^{(1)}>\rho$ and do not back-propagate it to stage-I. The loss function for stage-II is formulated as

To jointly optimize stage-I and stage-II, we define a loss function aiming at estimating parameters $\bm{\Theta}^{(1)}$, $\bm{\Theta}^{(2)}$, $\mathbf{W}$, and $\mathbf{b}$ by optimizing

where $h^{(1)}$ and $h^{(2)}$ are the fusion weights.

2 Reverse Connections

FCNs have shown good segmentation results in recent studies, especially for single organ segmentation. However, for multi-organ segmentation, lower layers typically lack semantic information, which may lead to inaccurate segmentation particularly for smaller structures. Therefore, we propose reverse connections which feed coarse-scale (high) layer information backward to fine-scale (low) layer for semantic segmentation of multi-scale structures, inspired by . This enables us to connect abstract high-level semantic information to the more detailed lower layers so that all the target organs have similar levels of details and abstract information at the same layer. The reverse connections framework for stage-I is shown in Fig. 3. Fig. 4 illustrates a reverse connection block. Let $\mathbf{R}_{n}$ denote the reverse connection map of the $n$-th convolutional layer in the backbone network, i.e. FCN in this study, where $\mathbf{C}_{n}$ is the output of the $n$-th convolutional layer. A convolutional layer (with $512$ channels by $3\times 3$ kernels) is added after $\mathbf{C}_{n}$, and a deconvolutional layer (with $512$ channels by $4\times 4$ kernels) is applied after $\mathbf{R}_{n+1}$. $\mathbf{R}_{n}$ is then obtained via an element-wise summation of these two maps. $\mathbf{R}_{7}$ is the output of a convolutional layer (with $512$ channels by $2\times 2$ kernels) grafted onto $\mathbf{C}_{7}$. Let $\mathbf{w}^{n}$ denote the corresponding weights for obtaining $\mathbf{R}_{n}$. Following , we add reverse connections from $\mathbf{C}_{4}$ to $\mathbf{C}_{7}$.

With these learnable reverse connections, the semantic information of the lower layers can be enriched. In order to drive learned reverse connection maps to produce segmentation results approaching the ground-truth, we make each reverse connection map associate with a classifier. As the side-output layers proposed in are designed for detection purposes, they are not suitable for our task. Instead we follow the side-outputs used in . More specifically, a convolutional layer (with $|\mathcal{L}|$ channels by $1\times 1$ kernels) is added on top of $\mathbf{R}_{n}$, whose output is denoted as $\mathbf{V}_{n}$, and followed by a deconvolutional layer (with $|\mathcal{L}|$ channels). We denote the weights of the $n$-th side-output layer by $\bm{\theta}^{n}$. The loss function for side-output layers $\mathcal{J}^{(s,1)}$ is defined as

In order to combine the learned reverse connection maps of fine layers and coarse layers, we add up the predictions (i.e., $\mathbf{V}_{n}$) of the reverse connection maps from high layer to low layer gradually. First, $\mathbf{V}_{6}$ is fused with a $2\times$ upsampling of $\mathbf{V}_{7}$ by an element-wisely addition. Then we follow the same strategy and gradually merge $\mathbf{V}_{5}$ and $\mathbf{V}_{4}$, as shown in Fig. 5. To obtain a fused activation map $\mathbf{A}^{(f,1)}=\{a_{i,l}^{(f,1)}\}_{i=1,...,H\times W,l=0,...,|\mathcal{L}|}$ from the activation map of both side-outputs (i.e., $\mathbf{A}^{(r,1)}$) and convolutional layers in the backbone network (i.e., $\mathbf{A}^{(b,1)}$), a scale function is adopted followed by an element-wise addition by

where $\mathbf{A}_{l}$ indicates the $l$-th channel of the activation map. $h^{(r,1)}_{l}$ and $h^{(b,1)}_{l}$ are fusion weights. Then the fused probability map, $\mathbf{P}^{(f,1)}=\{p_{i,l}^{(f,1)}\}_{i=1,...,H\times W,l=0,...,|\mathcal{L}|}$, can be obtained by $p^{(f,1)}_{i,l}=\sigma(a^{(f,1)}_{i,l})$. The final objective function for stage-I is defined by

where $h^{(b,1)}$, $h^{(s,1)}$ and $h^{(f,1)}$ are fusion weights, and

Note that in our full system with the two-stage organ-attention network and reverse connections, all the parameters are optimized simultaneously by standard back-propagation

3 Testing Phase

In the testing stage, given a slice $\mathbf{I}$, we obtain the stage-I and stage-II probability map by

where $f(\cdot,\cdot)$ is the network functions defined in Sec. 2.1. A fused probability map of $\mathbf{P}^{(1)}$ and $\mathbf{P}^{(2)}$ is then given by

The final label map $\mathbf{S}=\{s_{i}\}_{i=1,...,H\times W}$ is determined by $s_{i}=\arg\min_{l\in\mathcal{L}}p_{i,l}$.

Statistical Label Fusion Based on Local Structural Similarity

As described in Sec. 1, our OAN-RC is based on 2D images which is an extreme case of 3D anisotropic patches. In this section, we propose to fuse anisotropic information obtained from different viewing directions using isotropic 3D local patches to estimate the final segmentation. Let us denote the segmentation results by $\mathbf{S}^{j},(j=1,\ldots,M=3)$, which are obtained as described in Sec. 2.3 from the axial (Z), sagittal (X), and coronal (Y) OAN-RCs. Depending on the viewing directions, sectional images contain different structures and may have different texture patterns in the same organs. These differences can cause nonidentical segmentations by the deep network as shown in Fig. 6 in 3D. In addition, there is no guarantee of connectivity between neighbor slices by independent use of slices for training and testing. Possible naïve approaches for determining the final segmentation in 3D from the OAN-RC results can be boolean operations such as union or intersection. Majority voting (MV) is another candidate for efficient fusion, however, theses approaches assume the same global weights of OAN-RC results. From the observations that the performance level of segmentation, e.g. sensitivity, can be different from viewing directions for each organ, we set the performance level to be an unknown variable when computing the probability of labeling. This concept is similar to the label fusion algorithms using expectation-maximization (EM) framework such as STAPLE (simultaneous truth and performance level estimation) and its extensions .

Let us denote the true label of the $V$ by $\mathbf{T}$, which is unknown, and the unknown performance level parameter of segmentation by $\theta$. The segmentations from the deep networks $\mathbf{S}=\left\{\mathbf{S}^{j}|j=1,...,M\right\}$ are observed values. Under this condition, the basic EM framework is performed by following two steps in an iterative manner: 1) to compute $Q^{0}(\theta|\theta^{(k)})=E_{\mathbf{T}}\left[\ln L(\theta|\mathbf{S},\mathbf{T})|\mathbf{S},\theta^{(k)}\right]$ which is the expected value of the log likelihood, $\ln L(\theta|\mathbf{S},\mathbf{T})=\ln P(\mathbf{S},\mathbf{T}|\theta)$, under the current estimate of the parameters $\theta^{(k)}$ at $k^{th}$ iteration, and 2) to find the parameter $\theta^{(k+1)}$ which maximizes $Q^{0}(\theta|\theta^{(k)})$.

By assuming independence between $\mathbf{T}$ and $\theta$ in our problem, the second term $\sum_{\mathbf{T}}\ln P(\mathbf{T})P(\mathbf{T}|\mathbf{S},\theta^{(k)})$ in (13) becomes free of $\theta$ and the maximization step can be written as

Therefore, we redefine $Q^{0}(\theta|\theta^{(k)})$ as $Q(\theta|\theta^{(k)})=E_{\mathbf{T}}\left[\ln P(\mathbf{S}|\mathbf{T},\theta)|\mathbf{S},\theta^{(k)}\right]$.

The performance level parameter in this framework is a global property representing the overall confidence of deep network segmentation for the whole volume. However, it can also vary according to the voxel spatial locations via the local and neighbor structures as we use 2D slices for the initial segmentation. Therefore, we propose to combine local structural similarity shown from a specific viewing direction to the original 3D volume and the global performance level, conceptually similar to local weighted voting . We compute the probability of correspondence between 2D images and the 3D volume by structural similarity (SSIM) by

Considering the local image properties, the expectation of log likelihood function in our problem becomes

The global underlying performance level parameters of the deep network segmentations is defined as

where $\theta_{js^{\prime}s}$ is the probability of the voxel labeled as $s^{\prime}$ from the $j^{th}$ deep network with the current estimated performance value $\theta_{js^{\prime}s}^{(k)}$, when the true label is $s$.

To make the problem simple, we assume conditional independence between labeling and the original volume intensities. The labeling probability with the target image intensity then becomes

In the expectation step (E-step), we estimate the probability of voxelwise labels. Let us denote the probability that the true label of $i^{th}$ voxel is $s\in\mathcal{L}$ at the $k^{th}$ iteration by $\omega_{si}^{(k)}$. When the deep network segmentations $\mathbf{S}$ and performance level parameters at the $k^{th}$ iteration $\mathbf{\theta}^{(k)}$ are given, $\omega_{si}^{(k)}$ can be then described as

2 M-step

In the maximization step (M-step), the goal is to find the performance parameters, $\mathbf{\theta}$, which maximize (16) with the current given parameters. Considering each $\mathbf{S}^{j}$ and $\theta_{j}$ independently, the expectation of log likelihood function in (16) can be expressed with the estimated voxelwise probability in E-step. Then the performance parameter of each segmentation can be formulated to find the solution which maximizes the summation of voxelwise probability as

From the definition of $\theta$ in (17), the summation of probability mass function, $\sum_{s^{\prime}}\theta_{js^{\prime}s}^{(k)}$, must be $1$, and (22) becomes a constrained optimization problem which can be solved by introducing a Lagrange multiplier, $\lambda$. We then obtain the optimal solution by making the first gradient zero as

By applying the derivation of $Q$ in (16), (22) and (23), (24) becomes

By substituting the constraint of $\sum_{s^{\prime}}\theta_{js^{\prime}s}^{(k)}=1$, we can obtain the final optimal solution as

The two steps, (21) and (26), are then computed alternatively in the EM iterations until they converge. From the final values of (21), the final segmentation can be computed by graph-based approaches such as .

3 Parallel computing using GPUs

The fusion step can be efficiently computed in a parallel way on a GPU. The local structural similarity $\alpha_{i}^{j}$ of $i$-th voxel in $j$th deep network and priori $P(T_{i})$ can be computed for each voxel and saved as a pre-processing step. In the EM iterations, as shown in (21), the probability can be computed and updated for each structure at each voxel. In our implementation, a GPU thread is logically allocated for each voxel. However, to reduce the used memory and computation cost, the target volume of interest (VOI) for each structure $s$ is computed in an extended region as $\delta=4$ voxels for each direction from $V(\bigcup_{j}\mathbf{S}^{j}=s)$ in our implementation. For parallel computing, one CPU thread is allocated to a structure and launches a kernel of one GPU to compute EM iteration for each structure.

Experimental Results

We evaluated our methods on $236$ abdominal CT images of normal cases under an IRB (Institutional Review Board) approved protocol in Johns Hopkins Hospital as a part of the FELIX project for pancreatic cancer research . CT images were obtained by Siemens Healthineers (Erlangen,Germany) SOMATOM Sensation and Definition CT scanners. CT scans are composed of $(319-1051)$ slices of $(512\times 512)$ images, and have voxel spatial resolution of $\left([0.523-0.977]\times[0.523-0.977]\times 0.5\right)mm^{3}$. All CT scans are contrast enhanced images and obtained in the portal venous phase.

A total of $13$ structures for each case were segmented by four human annotators/raters, one case by one person, and confirmed by an independent senior expert. The structures include the aorta, colon, duodenum, gallbladder, interior vena cava (IVC), kidney (left, right), liver, pancreas, small bowel, spleen, stomach, and large veins. Vascular structures were segmented only outside of the organs in order to make the structures exclusive to each other (i.e. no overlaps).

As explained in Sec. 2, we used OAN-RCs for multi-organ segmentation whose backbone FCNs had been pre-trained by $PascalVOC$ dataset . From the possible variants of FCNs (e.g., FCN-32s, FCN-16s, and FCN-8s), which depend on how they combine the fine detailed predictions , we selected FCN-8s in this study because it captures very fine details in the $3^{rd}$ and $4^{th}$ pooling layer, and keeps high-level semantic contextual information from the final layer. Our algorithm was implemented and tested on a workstation with Intel i7-6850K CPU, NVidia TITAN X (PASCAL) GPU. With $236$ cases, the initial segmentations using OAN-RCs were tested by four-fold cross-validation. All the input images of OAN-RCs are $1.5$ times enlarged by upsampling, which lead to improved performance in our experiments.

In the fusion step, the average probability of $\mathbf{S}^{X},\mathbf{S}^{Y},\mathbf{S}^{Z}$ are taken as a priors in (21) and the initial performance levels $\theta_{js^{\prime}s}^{(0)}$ were computed by randomly selecting 5 cases and by comparing them to the ground-truth. To compute the local patch-based structural similarity in (15), patches of $(4.5\times 4.5\times 4.5)mm^{3}$ size cubes were used for 3D volume. Since CT voxels are not always isotropic and spatial resolutions can be different between scan volumes, we re-sampled the 3D patch with $0.5mm$ length cubic voxels so that the same size of $(9\times 9\times 9)$ 3D patches and $(9\times 9)$ 2D patches from all directions can be used for all cases in our experiments.

The final segmentation results using OAN-RC with local structural similarity-based statistical fusion (LSSF) were compared with the 3D-patch based state-of-the-art approaches, 3D Unet and hierarchical 3D FCN (HFCN) as well as 2D-based FCN, OAN and OAN-RC with majority voting (MV). For a quantitative comparison, we computed the well-known Dice-Sørensen similarity coefficient (DSC) and the surface distances based on the manual annotations as ground-truth. For a structure $s$, DSC is computed as ${{2V(\mathbf{S}=s\bigcap\mathbf{T}=s)}\over{V(\mathbf{S}=s)+V(\mathbf{T}=s)}}$ where $\mathbf{S}$ is the estimated segmentation and $\mathbf{T}$ is the ground-truth, i.e. manual annotations in this study. The surface distance was computed from each vertex of the ground-truth and to the estimates of our algorithms. Fig. 8 shows comparison results by box plots, while Tables 1 and 2 represent the mean and standard deviations for all the $236$ cases.

As shown in Fig. 8, the basic OAN-RC outperforms other state-of-the-art approaches and our local structural similarity-based fusion improves the results even more. We note that although DSC shows the relative overall volume similarity, it does not quantify the boundary smoothness or the boundary noise of the results. But evaluating the surface distances, see below, shows that our method works effectively for both the whole volumes and the boundaries of the organs.

Tables 1 and 2 represent the mean and standard deviations of performance measures for 13 critical organs. Similar to the box plots, they show that our OAN-RCs with statistical fusion improves the overall mean performance and also reduces the standard deviations significantly.

The OAN-RC training and testing can be computed in parallel for each view direction. In our experiments, the training took $40$ hours for $120,000$ iterations for $177$ training cases and the average testing time for each volume was $76.73$ seconds. The fusion time depended on the volume of the target structure, and the average computation time for $13$ organs was 6.87 seconds.

Discussion

Multi-organ segmentation using OAN-RCs alone, without the statistical fusion, gave similar or better performance compared with the state-of-the-art approaches summarized in . In the specific case of the pancreas, state-of-the-art methods showed (mean $\pm$ standard deviations) segmentation accuracies as $74.4\pm 20.2(\%)$ on 140 cases , $78.5\pm 14.0(\%)$ on 150 cases , $78.0\pm 8.2(\%)$ on $82$ cases and $75.74\pm 10.47(\%)$ (on the whole slice) versus $82.4\pm 5.7(\%)$ (reduced region of interest) on $82$ cases in terms of DSC. We cannot make a direct comparison because in these datasets CT images and manual segmentations (i.e. annotation) for the ground-truth are different from each other. But our OAN-RCs segmentations on our larger dataset shows similar or better performances in terms of DSC. Among target organs, our performance on structures such as gallbladder and pancreas, whose sizes are relatively small and have particularly weak boundaries improves significantly from using basic FCNs or using OANs without reverse connections.

Moreover, as shown in Sec. 4, our statistical fusion based on local structural similarity improves the overall segmentation accuracies in terms of both DSC and average surface distances. In particular, there are significant performance improvements for the minimum values as shown in Fig. 8, which helps explain the robustness of the algorithm. The differences can be depicted more clearly by visualizing the 3D surfaces as shown in Figs 10 - 11. The noise of the deep network segmentations is distributed over large regions, without much connectivity, and occasionally they show significantly different patterns. But our fusion step exploits structural similarity which outputs clean and smooth boundaries by effectively combining different information based on the local structure of the original 3D volume.

When applying our proposed method and interpreting the evaluation results, we must address several considerations.

As shown in our experiments, our proposed algorithm also outperforms 3D patch based approaches. But 3D (isotropic) patch-based approaches have several issues which make it hard to apply to this problem. To make bigger patch size, they require more parameters and hence require more training data or, if this is not available, significant data augmentation (.e.g, by scaling, rotation, and elastic deformation). In addition, there can be practical memory limitation on GPUs which restricts the expandable patch size. The limited patch size means that the deep networks receptive field sizes contains only limited local information which is problematic for multi-organ segmentation and the discontinuities between the patches also raises problems. It is possible that solutions to these three problems may make 3D patch based methods work better in the future. Unlike 3D approaches, the local structure-similarity used in our fusion method effectively combine the information from anisotropic patches to 3D at each voxel. Fig. 9 shows an example generated by our proposed algorithm, which is visually indistinguishable from manual segmentation for almost all target structures.

The ground-truth used in this study for training and evaluation was specified using manual annotations by human observers. It is well known that there can be significant inter-/intra-observer variations in manual segmentation. But, as explained before, the ground-truth was created by four human observers and checked by experts in a visual way, and we randomly divided testing groups in our 4-fold cross-validation to avoid biased comparison. However, it is still possible that inaccuracies due to human variability may affect the evaluation as well as the training. This can be further intensively explored as separate experiments.

Another possible consideration when applying the proposed approach is the image quality which can affect both of manual annotations and deep network segmentation results. Various factors such as spatial resolution, level of artifacts and reconstruction kernels should be considered. The dataset used in this study has been collected between $2005$ to $2009$ in the same institute with control over the scanning parameters. As explained in Sec. 4, the CT protocol is the portal venous phase and the spatial resolution is almost isotropic. But different scanning parameters and artifacts may affect our algorithms performance when applied to other datasets.

The same issues about manual segmentations and image qualities can be raised in general segmentation and evaluations. Specifically for our proposed approach, especially in the fusion step, the way of computing priori, $P(T)$, used in (21) can in practice affect the final segmentation. But considering that the deep network segmentation results from different viewing-directions are independently obtained, the mean can be accepted in general. However, if the deep network segmentations show clear tendencies towards over-estimation or under-estimation, then different types of models for priors may need to be used in order to improve the final result for practical applications.

One of the main advantages of our algorithm is the efficient computation time. The segmentation of $13$ organs of the whole volume takes similar to or less than 1 minute with better performance reported than the state-of-the-art methods . Hence our approach can be practically useful in clinical environments.

Conclusion

In this paper, we proposed a novel framework for multi-organ segmentation using OAN-RCs with statistical fusion exploiting structural similarity. Our two-stage organ-attention network reduces uncertainties at weak boundaries, focuses attention on organ regions with simple context, and adjusts FCN error by training the combination of original images and OAMs. Reverse connections deliver abstract level semantic information to lower layers so that hidden layers can be assisted to contain more semantic information and give good results even for small organs. The results are improved by the statistical fusion, based on local structural similarity, which smooths our noise and removes biases leading to better overall segmentation performance in terms of DSC and surface distances. We showed that our performance is better than previous state of the art algorithms. Our framework is not specific to any particular body region, but gives high quality and robust results for abdominal CTs, which are typically challenging regions due to their low contrast, large intra-/inter-variations, and different scales. In addition, the efficient computational time of our algorithm makes our approach practical for clinical environments such as CAD, CAS or RT.

References