A Hierarchical Algorithm for Multiphase Texture Image Segmentation

Abstract

Image segmentation is a fundamental task for many computer vision and image processing applications. There exist many useful and reliable models for two-phase segmentation. However, the multiphase segmentation is a more challenging problem than two phase segmentation, mainly due to strong dependence on initialization of solutions. In this paper we propose a reliable hierarchical algorithm for multiphase texture image segmentation by making full use of two-phase texture models in a fuzzy membership framework. Application of the new algorithm to the synthetic and real medical imaging data demonstrate more satisfactory results than existing algorithms.

1. Introduction

Image segmentation is a fundamental problem in image processing which is a prerequisite to high-level computer vision applications. It aims to divide an image representing a real scene or a synthetic one into classes or categories, corresponding to different objects and the background in the image. In the end, each pixel should belong to one class and only one. In other words, we look for a partition of the image into distinct segments, and each of them shares some features in common such as intensities, color, or texture. In particular, image texture defined by repeated patterns of intensities adds much complication in image processing tasks. A textured image often has several regions with different textures in existence, and the task of segmentation is to locate the texture boundaries, with which this paper is concerned.

Over the past two decades, a variety of different techniques have been developed to solve the problem of image segmentation, ranging from region growing and emerging [1], watershed algorithms [2], minimum description length criteria [3], and active contour models [4, 5] to Mumford-Shah energy minimization model [6].

Historically, image segmentation is known as the process to segment an image into two categories of regions: the foreground and the background. This process nowadays is referred to as a two-phase modeling, whereas multiphase modeling is specifically to deal with segmentation of more than two regions. Earlier work on segmentation attempts to detect the feature boundaries directly by edge detection [7–9]. These methods are susceptible to noise which is often present in real applications, and as such they are not suitable for processing either textured or noisy images, unless one applies it to the transformed image after applying a Gabor-type filter to the original one [10].

For nontexture images, the most influential model, known as MS-model, is proposed by Mumford and Shah [6], where the required boundary set Γ of features as well as the segmented image 𝑢 is determined from minimizing an energy functional. Unfortunately this elegant model is numerically difficult to realize, as such considerable effort has been made in order to alleviate this problem. For instance, Ambrosio and Tortorelli [11] have proposed a more solvable model by approximating the MS-model with elliptic functionals defined on Sobolev spaces. However, the most well-known paper based on [6] is the algorithm proposed by Chan and Vese [12]. Known as the CV-model, this was initially designed for two-phase segmentation of images of approximately piecewise constant intensities, using the framework of level set functions [13].

Segmentation of texture image is intrinsically more challenging than intensity-based ones. For texture images, detecting edges directly is problematic as it would treat texture patterns as image features. One of the strategies is to derive texture features from the image by first applying certain filters to it, and subsequently these features will be treated as a multichannel (i.e., color) segmentation problem. For instance, Sandberg et al. [14] have used Gabor filters for texture image segmentation by extending the original CV-model [12], and a similar strategy has also been adopted by a later work [15]. On the other hand, there have been other powerful models proposed to tackle texture segmentation of texture image directly [16–18]. Especially, both the region competition model of Zhu and Yuille [16] and the nonparametric model of Kim et al. [17] minimized mutual information as fidelity measure plus the regularization term. More recently Ni et al. [18] proposed a more reliable model based on Wasserstein distance as fidelity measure which leads to global optimal solutions to segmentation of texture image, where local histogram information is extensively used in order to divide the image domain into two regions in each of which the difference of the cumulative distribution function from its median is minimal.

Mory and Ardon [19] used the concept of fuzzy region competition to unify and extend it to model texture segmentation problems. Here, each pixel is assigned a probability, instead of a precise membership integer, of belonging to a particular region. The fuzzy region competition has been further generalized to deal with multiphase segmentation problems by Li and coworkers [20, 21], and in particular the multiphase texture segmentation model [20], referred to as LN-model, will be discussed in the following sections.

For the past decade, research in multiphase nontexture segmentation has been active with a focus on grayscale images. The work on multiphase texture segmentation is more recent.

A number of multiphase models have been proposed most of which are generalized from two-phase ones; see [22–30]. Gao and Bui [26] proposed a hierarchical model for multiphase segmentation problem for piecewise smoothing image. Jung et al. [31] formulated a multiphase segmentation model built upon the celebrated phase transition work of Modica and Mortola in material sciences. An efficient algorithm for minimizing the piecewise constant Mumford-Shah functional of image segmentation was proposed [32] based on the threshold dynamics of Merriman et al. [33] for evolving an interface by its mean curvature. More recent work include models based on shape and topological sensitivity [28] and 𝐻₁ regularization model [34]. Salah et al. [35] recently proposed to use a kernel function to map the original image into data of a higher dimension so that the piecewise constant model becomes applicable. Inclusion of shape constraints into the multiphase segmentation was also explored by Cremers [36]. Several models for multiphase segmentation of image frames have also been proposed [29, 37]. We also remarked that there is a very interesting unsupervised multiphase segmentation model from [30], which can automatically determine the number of regions during the segmentation process.

It has been observed that most generalized models work well only for a small number of problems, and robustness is a major issue. For instance, the above-mentioned CV-model [12] was generalized to multiphases in [22] and then refined in [23]. However as noted by [38, 39], the new model [23] has a strict requirement on its initial guess which contradicts the idea of automatic segmentation. Moreover, the phase number has to be 2^𝑛, where 𝑛 is the number of level set functions. To improve this model, the idea of hierarchical implementation of the robust CV-model [12] was first considered in [38] (with time-marching solver a.k.a. additive operator splitting (AOS)) and was further improved by Badshah and Chen [39] with a robust multigrid solver. Along this line, Ni et al. [40] employed the fast time-marching dual algorithm proposed by [41]. These improved hierarchical models have produced convincing results in segmentation of intensity images but not yet applied to segmentation of texture images.

In clear contrast to two-phase segmentation, multiphase texture segmentation is relatively understudied so far. Several feature-based models have been proposed to tackle the problem by working on features derived by different strategies. For instance, Aujol et al. [42] have proposed a model based on features derived from applying wavelet transforms to the image of interest. Similarly, Wei and Xin [43] proposed a supervised model based on contourlet features for aerial image segmentation. These feature-based models have limitation in choosing appropriate descriptors of texture. The most recent Li-Ng (LN) model [20] was based on mutual information which implies that there is no need to detect features first as is required by feature-based models. However, this model is not globally convex, as such the performance is sensitive to initialization; Li and Ng [20] have remarked that this occurs in segmentation problems with three or more phases. As such more advanced segmentation models are needed.

We finally remark that there have been recent advances in efficient solvers for variational models. Different from the past solvers such as gradient descent level set methods, a diversity of fast solvers have been proposed; these include graph cut [44], multigrid [39, 45], dual projection algorithm [41], as well as genetic algorithms [46]. Amongst these, alternating optimization strategies via the elegant dual projection model have been increasingly employed in the literature [20, 40].

The rest of the paper is organized as follows. Section 2 reviews some existing segmentation models, paying special attention to fuzzy region competition models which assign a probability value at each pixel rather than an integer for phases. Section 3 introduces our new hierarchical algorithm based on general two-phase models. Section 4 shows a series of experiments for comparisons and verifications. Some conclusions are drawn in Section 5.

2. Review of Some Segmentation Models

We shall review five models that are reliable for two-phase segmentation. The first two models are designed for non-texture segmentation, and the last three models are for texture segmentation; in particular the fourth model by Ni et al. is the most reliable as well as the most expensive.

Let Ω be a bounded open subset of ℝ² with Γ=𝜕Ω its boundary, and let 𝐼∶Ω→ℝ be the given grayscale image. The aim of multiphase segmentation is to partition Ω into 𝑁 regions Ω_𝑖, 1≤𝑖≤𝑁, where Ω_𝑖⋂Ω_𝑗=0,𝑖≠𝑗 and ∑^𝑁𝑖=1Ω_𝑖=Ω.

Mathematically, each phase may be characterized by a distinct constant as defined in [5] and in two phases with 𝑁=2, one may simply try to find a piecewise constant solution 𝑢 taking two constant values in {0,1} (or two other constants {𝛼₁,𝛼₂} for membership identity). However it is often advantageous to look for a solution 𝑢 that takes values in [0,1] in a fuzzy membership framework [19].

2.1. The Mumford-Shah-Based CV Model

A general form of two-phase segmentation problem can be represented [19] as min_{Γ,𝛼1,𝛼2}𝐸Γ,𝛼₁,𝛼₂=_Γ𝑑𝑠+𝜆₁_Ω1+𝜆₂𝑟₁𝛼₁,𝑥𝑑𝑥_Ω2𝑟₂𝛼₂,𝑥𝑑𝑥,(1) where Γ is the (unknown) union of all boundaries of Ω_𝑖, and 𝜆₁,𝜆₂ are two positive and suitable multipliers; often we choose 𝜆₁=𝜆₂=𝜆. The two functions 𝑟₁(𝛼₁,𝑥), and 𝑟₂(𝛼₂,𝑥) are respective similarity measures of pixels (or region errors) in domains Ω₁, and Ω₂ with their representative values.

The well-known CV model [12] chooses 𝑟_𝑖(𝛼_𝑖,𝑥)=(𝛼_𝑖−𝐼(𝑥))² in (1). Here, 𝑢 takes values in {𝛼₁,𝛼₂}. The model can be adapted to tackle the multiphase segmentation problem [38–40] when more constants are defined.

The level set formulation (as a means to rewrite the above integrals over Ω₁,Ω₂ as over Ω) is commonly used to represent this model, and subsequently the resulting partial differential equation (PDE) is solved by a gradient descent time-marching method, as with the celebrated CV model [12]. However, this model (1) is not convex (neither is the CV model), as such the level set methods need a reasonable initialization in order to avoid local minima. This prerequisite and the slow convergence are inherent weakness of this PDE-based method.

2.2. A Global Convex Model

Recently, Bresson et al. [47] have tried to reformulate the above problem into a convex one so that the global minimum becomes easier to compute. The idea is to introduce constraints and to convert the two-phase model into a convex total variational (TV) model as follows: min_{0≤𝑢≤1,𝛼1,𝛼2}𝐸𝑢,𝛼₁,𝛼₂=+_Ω𝜆||∇𝑢_Ω(𝑥)||𝑑𝑥+𝜆_Ω𝑢(𝑥)𝑟₁𝛼₁,𝑥𝑑𝑥(1−𝑢(𝑥))𝑟₂𝛼₂,𝑥𝑑𝑥.(2) It should be remarked that Bae and Tai [48] have proposed a convexified model for four-phase segmentation.

Further, the Chambolle’s pioneering work [41] of fast dual projection can be used to provide an elegant efficient solver for this equation. More specifically, after an auxiliary variable 𝑣 is added to it, (2) will be formulated as min_{0≤𝑢≤1,𝛼1,𝛼2}𝐸𝑢,𝛼₁,𝛼₂=_Ω||∇𝑢(𝑥)||𝑑𝑥+12𝜃_Ω+𝜆_Ω(𝑢(𝑥)−𝑣(𝑥))²𝑑𝑥+𝜆_Ω𝑣(𝑥)𝑟₁𝛼₁,𝑥𝑑𝑥(1−𝑣(𝑥))𝑟₂𝛼₂,𝑥𝑑𝑥,(3) where the convex form here ∫_Ω(𝑢(𝑥)−𝑣(𝑥))²𝑑𝑥 is to force 𝑢 and 𝑣 to be close to each other (while serving the purpose of decoupling the nonlinearity), 𝜃>0 is a small parameter to penalize the error between 𝑢(𝑥) and 𝑣(𝑥). One hopes that, on convergence (with smaller and smaller 𝜃), 𝑢=𝑣 is obtained. Then the model (3) becomes min_𝑢_Ω||∇𝑢(𝑥)||𝑑𝑥+12𝜃_Ωmin_{0≤𝑣≤1,𝛼1,𝛼2}(𝑢(𝑥)−𝑣(𝑥))²𝑑𝑥,(4)12𝜃_Ω+𝜆_Ω(𝑢(𝑥)−𝑣(𝑥))²𝑑𝑥𝑣(𝑥)𝑟₁𝛼₁,𝑥𝑑𝑥+𝜆_Ω(1−𝑣(𝑥))𝑟₂𝛼₂,𝑥𝑑𝑥,(5) where the minimization problem (4) can be efficiently solved by a fast dual projection algorithm [41] and (5) is solved explicitly. The derived solution is 𝑢(𝑥)=𝑣(𝑥)−𝜃div𝑤(𝑥),(6) where 𝑤 can be solved by a fixed point method as follows: 𝑤^𝑛+1=𝑤^𝑛+𝑑𝑡∇(div𝑤^𝑛−𝑣/𝜃)1+𝑑𝑡||div𝑤^𝑛−𝑣/𝜃||,(7) where 𝑑𝑡≤1/8 is some suitable time step. The solution of (5) can be derived as 𝑣(𝑥)=maxmin𝑢(𝑥)−𝜃𝜆𝑟₁𝛼₁,𝑥−𝑟₂𝛼₂,𝑥,1,0.(8)

Thus, we have discussed a simple algorithm for 𝑣(𝑥) and 𝑢(𝑥). We note that similar use of an intermediate variable such as 𝑣 here can be seen in other works for example, [49].

2.3. The Fuzzy Region Competition Form of Two-Phase Segmentation

The fuzzy region competition model by Zhu and Yuille [16] or Kim et al. [17] for texture segmentation chooses 𝑟_𝑖𝛼_𝑖,𝑥=log𝑃_𝑖𝐼(𝑥)∣𝛼_𝑖,(9) where 𝑃_𝑖(𝐼∣𝛼_𝑖) denotes the probability density function (usually a Gaussian) for a conditional probability. Although Zhu and Yuille have used parametric representation [16], and Kim et al. used nonparametric representation [17], both of them have used standard PDE solution methods. So these models again are not globally convex; as such solvers may suffer from getting stuck at local minima. Of note, if the probability density 𝑃_𝑖 is a Gaussian distribution and the variance is known, this model will then reduce to the aforementioned CV-model.

2.4. The Local Histogram-Based Model

Rather than the mutual information-based strategies, Ni et al. [18] proposed a new strategy using local histograms (i.e., Wasserstein distance, instead of image intensities) to define a region in order to segment texture images. First of all, a Wasserstein distance with exponent  1 (measuring the distance of two histograms 𝑃₁, and 𝑃₂ with the range [0,𝐿]) is defined as 𝑊₁𝑃₁,𝑃₂=^𝐿0||𝐹₁(𝑦)−𝐹₂(𝑦)||𝑑𝑦,(10) where for each 𝑦, 𝐹_𝑗(𝑦) is the cumulative histogram for 𝑃_𝑗; usually 𝐿=255 for intensity images. At each pixel 𝑥, an intensity-based histogram 𝑃_𝑥 of a small local ball 𝐵_𝑥,𝜏 of radius 𝜏 centered on 𝑥 is defined by 𝑃_𝑥(𝑦)=||𝑧∈𝐵^𝑥,𝜏∩Ω∶𝐼(𝑧)=𝑦||𝐹||𝐵_𝑥,𝜏∩Ω_𝑥||,(𝑦)=||𝑧∈𝐵^𝑥,𝜏∩Ω∶𝐼(𝑧)≤𝑦||||𝐵_𝑥,𝜏∩Ω||,(11) for each prescribed 𝑦. Then in the framework of (1), Ni et al. [18] chooses 𝑟_𝑖𝛼_𝑖,𝑥=^𝐿0||𝛼_𝑖−𝐹_𝑥(𝑦)||𝑑𝑦.(12) That is to say, an image is segmented according to how accumulated local histograms of a group of pixels are close to a fixed value. We remark that this approach so far has only been used in two-phase texture image segmentation. Since both 𝜏 and 𝐿 are involved, the model is expensive to apply.

2.5. The LN Multiphase Texture Model

Derived from fuzzy region competition method, the multiphase segmentation problem for texture images can be formulated through the concept of mutual information as follows [20]: 𝐸(Γ)=^𝑁_𝑖=1𝜇2_𝜕Ω𝑖𝑑𝑠−^𝑁_𝑖=1_Ω𝑖log𝑃_𝑖𝐼,Ω_𝑖𝑑𝑥,(13) where 𝑃_𝑖(𝐼,Ω_𝑖) is the nonparametric probability density function (pdf) that is estimated by Parzen window by the intensity values of pixels in region Ω_𝑖. In the case of 𝑁=2, the energy in the second term will be similar to the one by Kim et al. [17], who firstly proposed it for texture segmentation.

This energy can be reformulated in a total variational framework as follows: 𝐸(𝑈,𝑃)=^𝑁_𝑖=1_Ω||∇𝑢_𝑖||𝑑𝑥−^𝑁_𝑖=1𝜆_Ω𝑢_𝑖log𝑃_𝑖𝐼,𝑢_𝑖𝑑𝑥,(14) where 𝑈=(𝑢₁,…,𝑢_𝑁) is the fuzzy membership vector, and 𝑃=(𝑃₁,…,𝑃_𝑁) is the corresponding pdf vector and is then approximated by 𝐸(𝑈,𝑉,𝑃)=^𝑁_𝑖=1_Ω||∇𝑣_𝑖||𝑑𝑥+12𝜃_Ω𝑣_𝑖−𝜆_Ω−𝑢_𝑖²𝑑𝑥𝑢_𝑖log𝑃_𝑖𝐼,𝑢_𝑖𝑑𝑥.(15) By relaxing 𝑢_𝑁=1−∑^{𝑁−1𝑖=1}𝑢_𝑖, the above energy is further approximated by 𝐸𝑈,𝑉,𝑃=^𝑁−1_𝑖=1_Ω||∇𝑣_𝑖||𝑑𝑥+12𝜃_Ω𝑢_𝑖−𝜆_Ω−𝑢_𝑖²𝑑𝑥𝑢_𝑖log𝑃^𝑁𝐼,𝑢_𝑁𝑃_𝑖𝐼,𝑢_𝑖𝑑𝑥(16) subject to 0≤𝑢_𝑖≤1, for 𝑖=1,…,𝑁−1, where 𝑉=(𝑣₁,…,𝑣_𝑁).

Then minimization of this energy 𝐸(𝑈,𝑉,𝑃) was solved following the fast dual projection model iteratively and in each step:(i)for fixed 𝑃, 𝑉, 𝑢_𝑖=max{min{𝑣_𝑖−𝜃𝜆log(𝑃_𝑁/𝑃_𝑖,1},0}, 𝑖=1,…,𝑁−1;(ii)for fixed 𝑃 and 𝑈, 𝑣_𝑖=𝑢_𝑖−𝜃div𝑤_𝑖(𝑥), 𝑖=1,𝑁−1. Here, 𝑤_𝑖 can be solved by a fixed-point method similar to the two-phase problem, by initializing 𝑤^0𝑖=0, 𝑤^𝑖𝑛+1=𝑤^𝑛𝑖+𝑑𝑡∇div𝑤^𝑛𝑖−𝑣/𝜃1+𝑑𝑡||div𝑤^𝑛𝑖−𝑣/𝜃||,(17)(iii)for fixed 𝑈, and 𝑉, update 𝑃 as in [20].

In summary, the above two-phase models are reliable (when 𝑁=2), while the multiphase texture models are not (when 𝑁>2). Our idea below is to make full use of the reliable two-phase models for multiphase texture segmentation.

3. A Hierarchical Algorithm for Texture Segmentation

Since most multiphase segmentation models solve a nonconvex minimization problem, the solution methods can be seriously affected by initialization or easily get stuck in a minimizer, when the number of phases is larger than two. It should be remarked that there has been some promising progress recently [48] that can reformulate models into convex ones which have a global minimizer but no local minimizers. To improve on [20] when 𝑁>2, we pursue an alternative approach in this paper.

Motivated by [39] for a hierarchical algorithm and [1] for seeded region growing models, for non-texture images, we now propose a hierarchical algorithm for segmenting a multiphase texture image.

Our idea is the following. Assume that each phase contains a seeded point indicating the location of a phase; practically one simply clicks on a few points within the image to define seeds. We recursively apply a two-phase model (3) to split a partitioned region into two further subregions and generates an ordered binary tree to represent the structure of the image, as is illustrated in Figure 1. The entire image domain is initially defined as 𝑄₀, representing the initialization, which can be arbitrary, we may conveniently choose the image itself (normalized to range from 0 to 1) as initialization for 𝑣. Then, our two-phase model (i.e., mutual information or local histogram) is to segment the given image 𝐼, 𝑄₀ representing the entire domain, into two phases (a domain 𝑄₁ and its complement 𝑄₂) using fast dual projection strategy. The partitioned regions 𝑄₁ and 𝑄₂ will be stored as the tree nodes. For each of 𝑄₁ and 𝑄₂, there is a decision here to make to determine if they are further segmented. This can be done by checking how many seeds a subregion contains, and further segmentation is only carried out using the two-phase model again if more than 1 seed point is contained. Figure 1(a) illustrated the work flow and at the end the original image was successfully segmented into five regions (i.e., 𝑄₃, 𝑄₄, 𝑄₅, 𝑄₇, and 𝑄₈). Of note, region term can either be represented by mutual information or Wasserstein distance, or even an alteration between them may be implemented. Other strategies such as Mory and Ardon’s strategy [50] as cost function might be useful, or mixed form during the iteration might be useful.

(a) Example 1 (Figure 5(a)) demonstrating segmented phases

(b) Example 2 (Figure 3(a)) demonstrating segmented phases

(c) Entropy map of example 1

(d) Intrinsic texture map of example 1

(e) Entropy map of example 2

(f) Intrinsic texture map of example 2

Figure 1 

Two examples to illustrate our hierarchical algorithm.

For multiphase segmentation problems, it is always assumed that the number of phases are known and more than two, so at least one 2-phase segmentation will be unconditionally performed. For instance, for a 3-phase problem as illustrated in Figure 1(b), for 𝑄₁ and 𝑄₂, at least one of them has to be further segmented, but the program has to follow a rule to decide which one to continue. On the other hand, for the 5-phase problem Figure 1(a), both 𝑄₁ and 𝑄₂ should be further segmented.

In summary, our algorithm can be stated in the following steps.(i) initialization: Initialize 𝑄₀ and set seed points interactively.(ii) two-phase segmentation by(a) initialization: initialize 𝑢;(b) iterations: update 𝑣 by (8), where 𝑟 will be replaced with the strategy adopted, that is, mutual information strategy or local histogram; update 𝑢 by (6); update region term defined either as 𝑃_𝑖 or local histogram distance metrics;(c) termination;(iii) for each detected phase, run the previous step if there are two or more seeds in it.(iv) termination strategy: no further segmentation is possible.

Finally, we show some results on designing alternative and possible indicators that may be used to decide on further segmentation without using seeds. For example, region size, change of uniformity described by standard deviation, or even the difference between the child phases or between them and parent one, are candidate choices for this purpose. None of the indicators we tested are suitable for texture segmentation. Our analysis established that conventional measures that worked well for non-texture segmentation problems cannot be readily applied. In searching for new strategies, a measure that can reflect the uniformity of texture would be particularly attractive. From the many different options we tested, due to the limit of space, here we only discuss two of them: entropy and a new texture feature based on semilocal image information [51]. Maps of both of them were given in Figures 1(c), 1(d), 1(e), and 1(f), respectively. Values of mean and standard deviation (STD) of each region are presented in Tables 1 and 2. For the entropy one, an entropy map was produced by computing entropy values of a small window centred at each pixel in the image. By doing this, we wish we can reuse a concept similar to the work on non-texture images, unfortunately, this does work as stated. For the 3-phase problem, values of mean entropy ± STD are 4.91±0.29 and 5.85±0.18 for phases 𝑄₂ and 𝑄₃, respectively. If using STD as criteria, 𝑄₂ should be further split because of its nonuniformity, but this is incorrect. The second measure fails as well for this example. Using these measures for the 5-phase problem are also problematic mainly because that they do not monotonically decrease when the number of different textures becomes smaller. For instance, STD values of 𝑄₂ are lower than 𝑄₆. Combinations of these have also been considered, but this becomes a tuning process for threshold values from image to image. All of them are applicable only to a specific problem and a generic one proved to be nonexistent due to the complexity of texture. This automation is certainly one of the future research directions.

Table 1 

Example one for demonstration of stopping criteria (success: ✓ and failure: ×).

Table 2 

Example two for demonstration of stopping criteria.

4. Experimental Results

We first tested our algorithms on synthetic texture images and real medical images including optical coherence tomography (OCT) images. Note, here for the testing purpose, we set 𝜃=0.1, 𝜆=0.2, and 𝑑𝑡=1/8. For simplicity, we denote by:(i)M1—the method of CV model [12];(ii)M2—the method of mutual information by Li and Ng [20];(iii)M3—the method of local histogram by Ni et al. [40];(iv)M4I—the method of our framework using mutual information as fidelity term;(v)M4II—the method of our framework using local histogram as fidelity term.

4.1. Numerical Tests and Comparisons

Below, we use 2 sets of experiments to show that direct multiphase segmentation models can be sensitive to initialization, while our new algorithm is robust.

4.1.1. Comparisons of Two-Phase Models

In this first set of tests, we applied the CV-model (M1) and two texture segmentation models (M2 and M3) to synthetic and real texture images, see Figure 2. All three models performed reasonably well on a noisy synthetic grayscale image, but for texture image, the CV-model (M1) was no longer able to detect the object with textures. On the other hand, both texture models (M2 and M3) are good for two phases, while it appears that the mutual information one (M2) performs the best.

(a) Original image

(b) Result with M1

(c) Result with M2

(d) Result with M3

(e) Original image

(f) Result with M1

(g) Result with M2

(h) Result with M3

(i) Original image

(j) Result with M1

(k) Result with M2

(l) Result with M3

(m) Original Image

(n) Result with M1

(o) Result with M2

(p) Result with M3

Figure 2 

Simulation results on two-phase texture segmentation problems.

4.1.2. Segmentation of Multiphase Texture Images

In this second set of tests, we show that the M2 is good for three phases but can also fail for these phases if the initial guess is not chosen accurately enough, see Figure 3. Other two models M1 and M3 were not tested as the CV-model (M1) cannot cope with texture as shown in the previous, example while the M3 is not yet developed for multiphase segmentation.

(a) Original image

(b) Random initialisation

(c) Results from random initialisation

(d) User initialisation

(e) Results from user initialisation

(f) Original image

(g) Random initialisation

(h) Results from random initialisation

(i) User initialisation

(j) Results from user initialisation

Figure 3 

Results of three-phase segmentation by M2. Note: only two phases were segmented for (a).

We also applied our new segmentation algorithms M4I and M4II to the same three phase images and a five-phase segmentation problem, respectively, and the results are pleasing as expected; see Figure 4 for the results of three-phase problems and Figure 5 for the five-phase one.

(a) Original image

(b) Seed points

(c) Result with M4I

(d) Results with M4II

(e) Original image

(f) Seed points

(g) Result with M4I

(h) Results with M4II

Figure 4 

Results of 3-phase segmentation by our proposed strategies (M4I and M4II).

(a) Original image

(b) Seed points

(c) Result with M4I

(d) Result with M4II

Figure 5 

Results of five phase segmentation by our proposed strategies.

4.2. A Medical Imaging Application Segmentation of OCTs

Optical coherence tomography (OCT) is a noninvasive imaging technology that can reveal cross-sectional information of the biological tissues, and it is a powerful tool assisting diagnosis and management of a wide range of eye diseases. All OCT images are extremely noisy and none of the existing denoising models work well for such OCT images. In this section, we show a useful application of our proposed segmentation strategies to some OCT images, when treated as texture images. A typical OCT image of the retina demonstrates multilayers of the retina tissues, see Figures 6(a) and 6(e). The acceptable results by our M4I and M4II, as one can see from Figure 6, demonstrate the potential capability of our strategies for segmenting multiple retinal layers in the OCT images.

(a) Original image

(b) Seed points

(c) Result with M4I

(d) Result with M4II

(e) Original Image

(f) Seed points

(g) Result with M4I

(h) Result with M4II

Figure 6 

Results of segmentation of OCT images by our proposed strategies.

5. Conclusions

Multiphase texture segmentation problem represents a significant challenge in image processing. There is a lack of robust and reliable models in the literature. In this paper, we proposed a reliable hierarchical multiphase texture segmentation strategy based on recursive use of the two phase approaches. We have unified recently proposing two-phase strategies of representing region errors, that is, mutual information-based and Wasserstein distance-based, into the same framework. Our application to the synthetic and real medical imaging data demonstrated satisfactory results. One of our future work is to continue investigating alternative strategies for further segmentation (subdividing). We shall also investigate how to use our new algorithm for OCT layer segmentation for clinical research use.