Open access

Constrained Compound MRF Model with Bi-Level Line Field for Color Image Segmentation

Written By

P. K. Nanda and Sucheta Panda

Published: 24 October 2012

DOI: 10.5772/50475

From the Edited Volume

Advances in Image Segmentation

Edited by Pei-Gee Peter Ho

Chapter metrics overview

1,843 Chapter Downloads

View Full Metrics

1. Introduction

Image segmentation is a basic early vision problem which serves as precursor to many high level vision problems. Color image segmentation provides more information while solving high level vision problems such as, object recognition, shape analysis etc. Therefore, the problem of color image segmentation has been addressed more vigorously for more than one decade. Different color models such as RGB, HSV, YIQ, Ohta (I1, I2, I3), CIE(XYZ, Luv, Lab) are used to represent different colors [5]. From the reported study, HSV and (I1, I2, I3) have been extensively used for color image segmentation. Ohta color space is a very good approximation of the Karhunen-Loeve transformation of the RGB, and is very suitable for many image processing applications [1]. Image Modeling plays a crucial role in image analysis. Stochastic models, particularly MRF models, have been successfully used as the image model for image restoration and segmentation [2], [3], [4]. MRF model has also been successfully used as the image model while addressing the problem of color image segmentation both in supervised and unsupervised framework. Kato et al [6] have proposed a MRF model based unsupervised scheme for color image segmentation. In Kato 's method, the model parameters have been estimated using Maximum Likelihood criterion and the only parameter identified by the user is the number of class. This algorithm could be validated using different color textures and real images. Another color texture unsupervised segmentation algorithm has been proposed by Deng et al [7] and the method has been retermed as JSEG method. Recently, an unsupervised image segmentation algorithm has been proposed by Guo et al [8] where K-means has been used to initialize the classification in the classification of numbers. Very recently Scarpa et al. [13] have proposed a multiscale texture model and a related algorithm for the unsupervised segmentation of color images. In this scheme, the feature vectors have been collected and based on the feature vector the textures are then recursively merged giving rise to larger and more complex textures. This algorithm could successfully be tested on real world natural and remote sensing images. The model parameters can be estimated in both supervised and unsupervised framework [6].

In this piece of work, a Constrained Compound MRF model based color image segmentation scheme is proposed in unsupervised framework. We have used Ohta (I1, I2, I3) color space to model the color images. In the proposed scheme, the Constrained Compound MRF model parameters and the image labels are estimated concurrently. Since the image label estimates and the estimates of model parameters are dependent on each other, obtaining global estimates of label as well as model parameters is very hard. Hence, we have proposed a recursive scheme for estimation of image labels and model parameters. The recursive scheme yields partial optimal solutions as opposed to optimal solutions. The MRF model parameter estimation problem is formulated in Maximum Conditional Pseudo Likelihood (MCPL) framework and the MCPL estimates are obtained using homotopy continuation bases algorithm. The MCPL estimation strategy results in a set of nonlinear equations which need to be solved to determine the model parameter estimates. Determination of the estimates is tantamount to determine the zeros of the unknown function. Homotopy continuation methods [14], [15] are globally convergent methods that have been used to trace the zeros of a function and hence determines the solution of functions. We have developed the fixed point based homotopy continuation method to estimate the model parameters. The image label estimation problem is formulated in Maximum a Posteriori (MAP) framework and the MAP estimates are obtained using the proposed hybrid algorithm [10]. The proposed supervised algorithm has been successfully tested on different images, however, for the sake of illustration we have presented three results and a comparison is made with [9].

Advertisement

2. MRF model

MRF theory is a branch of probability theory for analyzing the spatial or contextual dependencies of physical phenomena. It is used in visual labeling to establish probabilistic distributions of interacting labels.

2.1. Neighborhood system and cliques

The sites in S are related to one another via a neighborhood system. A neighborhood system for S is defined as

N={Ni/iS}E1

Where Ni is the set of sites neighboring i. The neighboring relationship has the following properties:

  1. A site is not neighboring to itself: i N i

  2. The neighboring relationship is mutual: i N i   N i  

Ni={iS/[dist((xi,yi),(xi,yi))]2r,iiE2

For a regular lattice S, the set of neighbors of i is defined as the set of sites within a radius of r from i.

Where dist (A, B) denotes the Euclidean distance between A and B and r takes an integer value. The Fig 1 shows (η1) the first order and second order neighborhood system (η2).

Figure 1.

(a) Figure showing first order (η1), second order (η2) and third order(η3)neighborhood structure (b) Cliques on a lattice of regular sites.

The pair (S, N) = G constitutes a graph in the usual sense; s contains the nodes and N determines the links between the nodes according to the neighboring relationship. A clique c for (s, N) is defined as a subset of sites c={i, i), or a triple of neighboring sites c = {i, i, i’’), and so on. The collections of single-site, pair-site and triple-site cliques will be denoted by C1, C2, and C3 respectively, where

C1={i/iS}E3
C2={{i,i'}/i'Ni,iS}E4
C3={{i,i',i''}/i,i',i''Sareneighborstooneanother}E5

The sites in a clique are ordered, and {i, i} is not the same clique as {i, i}, and so on. The collection of all cliques for (S, N) is

C=C1C2C3.....E6

The type of a clique for (S, N) of a regular lattice is detetrmined by its size, shape and orientation. Fig. 1 shows the clique types for the first order and second order neighborhood systems for a lattice [2] [3].

Let Z={Z1,Z2,...,Zm} be a family of random variables defined on the set S, in which each random variable Zi takes a value zi in L. The family Z is called a random field. We use the notion Zi=zi to denote the event that Zi takes the value zi and the notion (Z1=z1,Z2=z2,....,Zm=zm).

To denote the joint event. For simplicity a joint event is abbreviated as Z = z where z={z1,z2,...} is a configuration of z, corresponding to realization of a field. For a discrete label set L. the probability that random variable Zi takes the value zi is denoted P(Zi=zi), abbreviated P(zi) and the joint probability is denoted as P(Z=z)=P(Z1=z1,Z2=z2,...,Zm=zm) and abbreviated P(z).

F is said to be a Markov Random Field on S with respect to a neighborhood system N if and only if the following two conditions are satisfied:

P(Z=z)>0,zZ(Positivity)E7
P(zi/zSi)=P(zi/zNi)(Markovianity)E8

Where S-i is the set difference, zS-I denotes the set of labels at the sites in S-i and

zNi={zi'/i'Ni}E9

Stands for the set of labels at the sites neighboring i.

The positivity is assumed for some technical reasons and can usually be satisfied in practice. The Markovianity depicts the local characteristics of Z. In MRF, only neighboring labels have direct interactions with each other[2][3].

The concept of MRF is a generalization of that of Markov processes (MPs) which are widely used in sequence analyisis. An MP is defined on a domain of time rather than space. It is a sequence of random variables Z1, Z2, …., Zm defined in the time indices 1, 2, …, m. It is generalized into MRFs when the time indices are considered as spaial indices.

Advertisement

3. Compound Markov Random Field (COMRF) model

Capturing salient spatial properties of an image lead to the development of image models. MRF theory provides a convenient and consistent way to model context dependent entities for e.g. image pixels and correlated features [6]. Though the MRF model takes into account the local spatial interactions, it has its limitations in modeling natural scenes of distinct regions. In case of color models, it is known that there is a correlation among the color components of RGB model. In our formulation, we have decorrelated the color components and introduced an interaction process to improve the segmentation accuracy. We have employed inter-color-plane interaction (Ohta I1, I2, I3 color model) process which reinforces partial correlation among different color components.

In this work, a compound MRF model has been proposed and the proposed model is based on the following notion. The prior MRF model takes care of (i)Intra-color-plane I1 or I2 or I3 I1, I2, and I3 entities of each color plane(ii)Inter-color-plane interactions of pixels of different color planes for e.g. I1 and I2, and I3. The MRF prior model takes care of the spatial interactions in any given color plane and also interaction of a pixel of a given color plane with the pixel of other color planes. Thus the intra color plane and inter color plane interactions could be modeled by the compound MRF model. Motivation behind this modeling is as follows. It is known that strong correlation exists among different color planes of RGB model and therefore not suitable for image segmentation. On the other hand Ohta model is suitable for image segmentation because of the existing weak correlation among color planes. In order to develop an appropriate color model, we develop a model with controlled correlation among the different color planes. Therefore, the a prior compound MRF model takes care of the controlled correlation among the different planes of Ohta colorspace. The degree of correlation is controlled by the associated parameters of the clique potential function I1, I2, I3. The values of these parameters are quite low and hence provide a controlled weak correlation among the inter planes making it suitable for image segmentation.

We assume all images to be defined on discrete rectangular lattice MxN. In the following the Compound MRF model is developed. Let the observed image X be modeled as a random field and x is a realization which is the given image. Let Z denote the label process associated with the segmented image Fig.2(a) shows the three planes of Ohta color model. Each color plane is modeled by a MRF model. Let L denote the number of labels. For a given plane for example Z, if the spatial interactions are modeled by MRF, then the prior probability distribution P(Z) is Gibbs distributed and can be expressed as

P(Z1=z1|θ)=(1/z)eU(z1,θ)E10

where Z=zeU(z,θ) is the partition function, U(z1,θ) is the energy function and is of the form U(z',θ)=cCVc(z',θ) being referred to as cliqe potential function, θ denotes the clique parameter vector. Analogously the spatial interactions of I2 and I3 planes can be defined. This prior MRF model taking care of all the three spatial planes would result in the energy function of the following form

P(Z=z|θ)=(1/z)eU(z,θ)z=zeU(z,θ)U(z,θ)=cCVc(z,θ)E11

where, Vc(z,θ) denotes the clique potential function for the three spatial planes I1, I2 and I3 respectively. However, the model is not complete for the color model. We model Z as a compound MRF, where the spatial interactions of individual color planes are taken care together with the inter color plane interactions of pixels. The inter color plane interactions of pixels of one plane with the other is shown in Fig.1(a). For the sake of illustration, Fig.1(b) shows interaction of (i, j) {th} pixel of I2 plane with the pixels of I1 plane with the first order neighbourhood structure in the inter color plane direction. If this inter color plane interactions need to modeled with the MRF prior, we can express

P(ZI2i,j=zI2i,j/ZI1k,lzI1k,l,(k,l)I1)=P(ZI2i,j=zI2i,j/ZI1k,l=zI1k,l,(k,l)(i,j),(k,l)ηI1i,j)

Figure 2.

(a) I1, I2, I3 Plane Interaction (b) Interaction of one pixel of I1 –plane with I2 -plane.

Let z denote the labels of pixels taking care of all three color planes. In otherwords, z denotes the labels for pixels of the color image. For example, z{i, j} corresponds to the (i, j) th pixel label consisting of three color components. The prior probability of z has been contributed by the intra color plane interactions and inter color plane interactions of pixels. hence, the prior model of z consists of the clique potential functions Vcs(z) and Vct(z) corresponding to intra color plane interactions and inter color plane interactions respectively. The vertical and horizontal line fields for different color planes (k=1, 2, 3) are denoted as v{k} and h{k} respectively. The horizontal and vertical line fields are defined as follows. Let fv(zki,j,zki,j1) for the kth color plane be defined as fv(zki,j,zki,j1)=|zki,jzki1,j| fh(zki,j,zki,j1)>threshold. Vertical line field for each plane is set i.e.

v{i, j} {k}=1 for k=1, 2, 3, if fv(zki,j,zki,j1)>threshold, else v {i, j} {k}=0. Similarly, in case of horizontal line field let fh(zki1,j,zki,j) be defined as fh(zki,j,zki,j1)=|zki,jzki1,j|.

Horizontal line field for k th plane is set, i.e. h {i, j} {k} =1 for k=1, 2, 3, if else h {i, j} {k} =1. Since the compound MRF model takes care of intra color plane as well as inter color plane interactions the prior probability distribution is given by (10), where the energy function can be expressed as,

U(z,θ)=Us(z,θ)+Ut(z,θ)E12

Where,

Us(z,θ)=i,jVcs(z,θ)E13
Ut(z,θ)=i,jVct(z,θ)E14

Here, Us(z, θ) and Ut(z, θ) refers to the energy function of intra-color-plane and inter-color-plane respectively. Vcs (zi, j) corresponds to the intra-color-plane pixels and Vct (zi, j) corresponds to inter-color-plane pixels. Let hsk for k=1, 2, 3 denote the horizontal line field for each color plane in intra-color-plane and htk for k=1, 2, 3 denote the vertical line fields for inter-color-plane directions. Thus the compound MRF model will have the energy function given by (12). Equation (13) can be written as,

Vcs(zi,j)=k=13αk[(zki,jzki,j1)2(1vki,j)+(zki,jzki1,j)2(1hki,j)]+βk[vi,j+hi,j]E15

Here, z1, z2, z3 correspond to I1, I2, I3 planes respectively. The equation (14) can be written as,

Vct(zi,j)=α1[(z1i,jz2i,j1)2(1v1i,j)+(z1i,jz2i,j1)2(1h1i,j)]+β1[v1i,j+h1i,j]+α2[(z2i,jz3i,j1)2(1v2i,j)+(z2i,jz3i,j1)2(1h2i,j)]+β2[v2i,j+h2i,j]+α3[(z3i,jz1i,j1)2(1v3i,j)+(z3i,jz1i,j1)2(1h3i,j)]+β3[v3i,j+h3i,j]E16

Where z1 denotes the interaction between I1-I2 color planes, z2 denotes the interaction between I2-I3 color planes and z3 denotes the interaction between I3-I1 color planes respectively. Here we have assumed, α1=α2=α3=α and β1=β2=β3=β. The [α,β]T is the set of unknown parameter vector that are selected on ad hoc basis. Since the line fields correspond to the edge pixels and in turn the boundary of a given segment. The similarity measure in case of boundary pixels for k=1, 2, 3 is not required and hence for boundary pixels, i.e. when hi, j=1 and vi, j=1, for k=1, 2, 3 the clique potential function of (16) consists of only the penalty function. Therefore the boundary pixels should not participate in the formation of regions with similarity measure.

Advertisement

4. Constrained Markov Random Field (MRF) model

In probability theory, a martingale is a stochastic process (i.e., a sequence of random variables) such that the conditional expected value of an observation at some time t, given all the observations up to some earlier time s, is equal to the observation at that earlier time s. Precise definitions are given below.

Originally, martingale referred to a class of betting strategies popular during 18th century, in France. The simplest of these strategies was designed for a game in which the gambler wins his stake if a coin comes up heads and loses it if the coin comes up tails. The strategy had the gambler double his bet after every loss, so that the first win would recover all previous losses plus win a profit equal to the original stake. Since as a gambler's wealth and available time jointly approach infinity his probability of eventually flipping heads approaches 1. The martingale betting strategy was seen as a sure thing by those who practiced it. Of course in reality the exponential growth of the bets would eventually bankrupt those foolish enough to use the martingale for a long time. The concept of martingale in probability theory was introduced by Paul Pierre Lévy, and much of the original development of the theory was done by Joseph Leo Doob. Part of the motivation for that work was to show the impossibility of successful betting strategies.

A discrete-time martingale is a discrete-time stochastic process (i.e., a sequence of random variables X1, X2, X3 that satisfies for all n,

E(|Xn|)<E(Xn+1/X1,X2,X3,...,Xn)=Xn

i.e., the conditional expected value of the next observation, given all of the past bservations, is equal to the last observation.

Capturing the salient spatial properties of an image lead to the development of image models [3]. Though the MRF model takes into account the local spatial interactions, it has its limitations in modeling natural scenes of distinct regions. In order to incorporate a stronger local dependence, we constrain this model based on the notion of martingale.The motivation behind the new model is as follows.

MRF model takes care of the local spatial interactions, nevertheless it has limitation in modeling natural scenes. In the following we propose new model with a view to take care of intra as well as inter plane interactions. In this research work, we employed the notion of martingale to reinforce the local dependence. Let Z(i),i=1,2,......n be a martingale sequence, namely for all i=1,2,......n E[|Z(n)|]<

and E[Z(n+1)/Z(1),.....Z(n)]=Z(n). Now, let Z1,Z2,.....Zn be the random variables associated with the image of size n=N2 and G is the predefined number of class labels. Therefore, E[Zi,j/Zk,l,k,li,j]=Zi1,j for any k,lηi,j, where ηi,j is the neighborhood of i,j. Consider,

E[Zi,j|Zk,l,k,li,j]=zi,jLzi,jP[Zi,j=zi,j|Zk,l=zk,l,k,li,j]E17

Assuming further that Z is a Markov process, we have

E[Zi,j|Zk,l,k,li,j]=zi,jLzi,jP[Zi,j=zi,j/Zk,l=zk,l,k,lηi,j]zi,jLzi,jP(Z=z)zi,jLP(Z=z)E18

Since Z is a MRF,

E[Zi,j|Zk,l,k,li,j]=zi,jLzi,jP(Z=z)zi,jLP(Z=z)E19

Since Zi, j is a martingle sequence E[Zi,j|Zk,l,k,li,j]=zk,lk,lηi,j

zk,l,k,lηi,j=zi,jLzi,jeU(z)zi,jLeU(z)E20

Considering first order neighbourhood and choosing one of the neighbourhood pixels for example zi-1, j, equation(20) can be expressed as

zi1,j=zi,jLzi,jeU(z)zi,jLeU(z)

Instead of taking a given pixel from the neighbourhood zi-1, j, we take the average of the neighborhood pixels. The a priori model of Z takes care of this constraint and the U(Z) is modified as (for∀ (i, j))

U(z)=i,jU(zi,j)+λc{zi,javgzi,jLzi,jeU(z)zi,jLeU(z)}2E21

Where zi,javg=zi,jLzi,jeU(z)zi,jLeU(z) and λc is the constrained model parameter. The energy function consists of two terms

U(z)=cCinVc(zs1,zs2,zs3)+cCirVc(zt1,zt2,zt3)E22

Where Vc(zs1,zs2,zs3) and Vc(zt1,zt2,zt3) are given by (15) and (16) respectively.

4.1. Constrained Compound Markov Random Field (CCMRF) model

The notions of the Constrained model has been fused with the notion of Compound Model to develop a new model known as Constrained Compound Model [10].

The model is given by

Usc(Z)=i,jU(zi,j)+λc{zi,javgzi,jLzi,jeU(zi,j)zi,jLeU(zi,j)}2E23

Where,

Usc(Z)=Vcs(zi,j)+λc{zi,javgzi,jLzi,jeU(zi,j)zi,jLeU(zi,j)}E24

Where Usc denote the energy function corresponding to intra color plane interactions and Vcs(zi, j) is defined by (15). Where, Zi,javg=zi,jLzi,jeU(zi,j)zi,jLeU(zi,j) and λc is the constrained model parameter. The energy function taking care of both intra-color-plane and inter-color-plane interactions with intra plane constraints is given by

U(Z)=Usc(zi,j,θ)+Utc(zi,j,θ)E25

Where Usc(zi,j,θ) is defined by (24) and Utc(zi,j,θ) is defined by (14).

Vcs(zi,j) and Vct(zi,j) are given by (15) and (16) respectively.

Advertisement

5. Unsupervised framework

In unsupervised scheme, the MAP estimates of the labels and the estimates of the model parameters are carried out concurrently. Thus, an estimation strategy need to be developed, which using the observed image, X, will yield an optimal pair (Zopt, θopt). The following joint optimality criterion is considered,

(zopt,θopt)=argz,θmaxP(Z=z|X=x,θ)E26

The estimated pair satisfying (26) is the global optima of P(Z=z/X=x, θ) with respect to Z and θ. Since both the entities Z and θ are unknown, and interdependent the problem is a very hard problem. Therefore, it is necessary to opt for strategies for suboptimal solution. In (26), z, θcould be viewed as a set of parameter of the given function P(Z=z/X=x, θ). For such kind of problems in deterministic framework, Wendell and Horter have proposed an alternate approach that would yield suboptimal solutions instead of optimal solution. Their approach is based on splitting the variables followed by recursively estimating the parameters. The final estimate in this process is called as the partial optimal solution. In our case, in stochastic framework, we in the same spirit venture to split the original problem into estimation of labels (z) and parameters estimate θ to obtain the partial optimal solutions. The splitting of the variables can be expressed as follows

(z*)=argzmaxP(Z=z|X=x,θ*)E27
(θ*)=argθmaxP(Z=z*|X=x.θ*)E28

These partial optimal solutions Z *and θ* are not global maxima, rather they are almost always local optimal solutions. But with θ=θ*, the estimate z* is global optimal satisfying equation (27) and analogously for z=z*, θ* is global optimal satisfying equation (28). Since neither θ* nor z* is known, a recursive scheme is adopted where the model parameter estimation and segmentation is alternated. Let at the kth iteration θk=[αk,φk]T be the estimate of model parameters and zk be the estimate of the labels of the observed image. We adopt the following recursion

(zk+1)=argzmaxP(Z=z|X=x.θk)E29
(θk+1)=argθmaxP(Z=zk+1|X=x,θ*)E30

The first problem of equation (29) is solved using Bayesian approach [2]. The optimal value of θ k is obtained by the proposed Homotopy Continuation method [6]. The MAP estimates are obtained by the proposed hybrid algorithm. One estimate of zk and θk constitute one combined iteration. this recursion is continued for finite number of steps to obtain zk and θk. Thus, the partial optimal solutions are obtained.

Advertisement

6. Image label estimation

The segmentation problem is cast as the pixel labeling problem. Each pixel can assume a label from the set of labels {0 − L}. In a given image of size L = M1 x M2, let Zi, j denote the random variable for (i, j)th pixel, (i, j) є L = M1 x M2. Z denotes the label process and z denotes a realization of the process. The label estimates z^ is obtained by maximizing the posterior probability P(Z=z|X=x,θ). Thus, the optimality criterion can be expressed as follows,

z^=argzmaxP(Z=z|X=x.θ^)E31

where, θ denotes the associated parameter vector of the double MRF model Z. Since z is unknown the above equation can not be computed. So, by using Baye’s theorem, hence (31) can be expre ssed as

z^=argzmaxP(X=x|Z=z,θ)P(Z=z)P(X=x|θ)E32

The observed image X is given and hence the denominator P(X=x|θ) of (32) is a constant quantity. P(Z = z) is the a priori probability distribution of the labels. The degradation process is assumed to be Gaussian and hence P(X=x|Z=z,θ) of (32) can be written as P(X=x|Z=z,θ)=P(X=z+w|Z,θ)=P(W=xz|Z,θ). Since, W is a Gaussian process, and there are three spectral components present in a color image, we have,

P(W=xz|Z,θ)=1(2π)ndet[K¯]12(xz)TK¯1(xz)E33

Where K is the covariance matrix. Hence, this minimization can be expressed as,

z^=argzmini,jk=13(x(i)z(i))22σ2+vcs(zki,j)+vct(zki,j)E34

Vcs(zi,j) and Vct(zi,j) are given by (15) and (16) respectively. Solving (34) yields the MAP estimates of the image labels and hence segmentation. The color image has three spectral components xk, zk, k=1, 2, 3, Vc is the clique potential function for all the three spectral components.

Advertisement

7. Model parameter estimation

We estimate the a priori model parameter using the ground truth image z. The associated MRF parameters of this ground truth image is θ. We also assume the number of labels associated with the original image to be known. The parameter estimation problem is formulated using Maximum Likelihood criterion. Here the image label available at the (k+1)th iteration is used to estimate θ at (k+1)th iteration. Therefore, the problem can be stated as the following

φk+1=argmaxP(Z=zk+1/θ)E35

Since, Z is a MRF, we have,

φk+1=argmaxθexp(U(zk+1,θ))ζexp(U(ζ,θ))E36

where ζ ranges over all realizations of the image z. Because of the denominator of (36), computation of the joint probability P(Z=zk+1/θ) is extremely difficult task. We maximize the pseudolikelihood function P^(Z=zk+1/θ) instead of the likelihood function P(Z=z/θ) where

(i,j)LP(Zi,j=zi,jk+1/Zm,n=zm,nk+1,(m,n)ηi,j,θ)=P(Z=zk+1/θ)E37

From the definition of marginal conditional probability, we can write

(i,j)LP(Zi,j=zi,jk+1/Zk,l=zk,lk+1,(k,l)(i,j),(i,j)L,θ)=P(Z=zk+1/θ)zi,jMP(X=x/θ)E38

Because of MRF assumption,

(i,j)LP(Zi,j=zi,jk+1/Zm,n=zm,nk+1,m,nηi,j,θ)=exp(cCVc(zk+1,θ))zi,jMcCVc(zk+1,θ)E39

Substituting equation (39) in (37) we have

P^(Z=zk+1/θ)exp(cCVc(zk+1,θ))zi,jMexp(Vc(zk+1,θ))E40

Therefore, the maximization problem (41) reduces to

argmaxθP^(Z=zk+1/θ)=argmax(i,j)Lexp(cCVc(zk+1,θ))zi,jMexp(cCVc(zk+1,θ))E41

In (41), the summation is over all possible labels M. (41) is highly nonlinear in nature and no a priori knowledge of the solution is available. Solving the resulting non-linear equations is hard and hence we developed a globally convergent based Homotopy Continuation method. We carry out the maximization process and obtain the estimate of parameter vector θ with the help of homotopy continuation method based algorithm.

7.1. Salient steps of the unsupervised algorithm

  1. Initialize parameter vector as θ0, pixel label estimates z0 for k=0, 1, 2, ..., N do

  2. Using θk, observed image x and initial segmented image zk, obtain the MAP estimate of the labels z^k+1

  3. With z^k+1, obtain the MCPL estimate of the parameter vector θ^k+1, using homotopy continuation based algorithm

  4. Compare θ^k+1with the previous estimate of θ^k, if |θ^k+1θ^k|<threshold, set θk+1=θkgo to step 2 else go to step 5

  5. Set estimate of parameter vector θ*=θ^k+1

  6. Estimate z* (segmented image) using θ*, z^k+1and observed image x

Advertisement

8. Parameter estimation using homotopy continuation method

8.1. Homotopy continuation method

Often, a wide variety of practical problems reduces to finding solution to a system of non-linear equations. The problem becomes difficult when we have little knowledge about the solutions of the system. In such situations, the popular Newton algorithm may fail to converge to a solution. Such examples can be found in [19]. Therefore, we need a method which, irrespective of the starting point always converges to a solution of the given system of equations. Homotopy continuation methods under some conditions always converges to a solution with probability one. Such methods are called globally convergent homotopy continuation methods [14]. The homotopy function is defined as follows :Let X, Y be two topological spaces and I be the unit interval λ/0λ1. The two maps f, g be maps from a space X to a space Y f,g:XY, then f is said to be homotopic to g if there exists a map H:XYsuch that H(x, 0) = f(x) and H(x, 1) = g(x) for x[0,1], such a map H is called a homotopy from f to g. In the above definition, H represents a continuous deformation of the map f to g as the parameter λ is varied from 0 to 1. There is no unique homotopy map that will continuously deform from a trivial map to any map. Depending upon the problem at hand the path has to be accurately tracked and hence, a suitable homotopy function has to be chosen for the existence of a path leading to the solution. The commonly used Homotopy maps are (i) Linear Homtopy (ii) Newton Homotopy (iii) Fixed Point Homotopy.

It is clear from Section 7 that the parameter estimation problem has been reduced to maximization of (41) with respect to θ. Towards this end let

f(θ)=θ{log[P^(X=xk+1/Y=y,θ)]}E42

Now the homotopy method is employed to solve f(θ)=0. In the following, we develop a general framework for solving f(θ)=0 using homotopy continuation method where θ is the unknown parameter vector to be determined.

In the continuation method we need to trace the homotopy path from a solution of a known system to that of the desired solution. In this regard, we have considered the fixed point homotopy map [14] which offers the advantage of arbitrary starting point for the path. This fixed point map is given by

h(θ,λ,q)=λf(θ)+(1λ)(θq)E43

where 0λ1 and q is an arbitrary starting point. Here the predictor-corrector method is employed to track the path defined by the homotopy in (43). The procedure can be briefly outlined as follows:

Let (θk,λk,θk1) be a point that satisfies (43). Therefore, the point considered is on the path. Tracking the path involves computing the adjacent point on the path. This is determined in the following way. Increment λk by some small value Δλ thus giving the next point λk+1=λk+Δλ and evaluate equation (43) at (θk,λk+1,θk1). If the value of the map h(θk,λk+1,θk1) is not equal to zero, then the point (θk,λk+1,θk1) is not on the path. Sinceh(θk,λk+1,θk1)0, we try to obtain an estimate of θk, say θ^k corresponding to λk+1 such that h(θk,λk+1,θk1)0. To achieve this one could use Newton's algorithm, namely,

θ^i+1k=θ^ikJθ^1[h(θ^ik,λk+1,θk1)]h(θ^ik,λk+1,θk1)E44

Where the superscript i denotes the ith Newton iteration and is the inverse of the Jacobian of h with respect to the coefficient of the parameter vector θ. But if θ^0k is too far from θ^k the value which makes, h(θ^k,λk+1,θk1)0 then (44) may not converge. To improve the convergence of (44), we select the initial point as θ^0k=θk. A further improvement in the convergence is obtained by considering

θ^0k=θ^kΔλJθ^1[h(θk,λk+1,θk1)]hλh(θk,λk+1,θk1)E45

The derivation of (45) is analogous to the derivation of Stonick and Alexander [15] for our homotopy map (43). Equation (45) corresponds to the prediction of the next point by taking a step in the direction of the path's slope. For the fixed point homotopy map considered, (45) becomes

θ^0k=θkΔλ{I1(λk+Δλ)I(1(λk+Δλ))2[Fθ1(θk)(λk+Δλ)+I1(λk+Δλ)]1}{f(θk)(θkθk1)}E46

Where I is the identity matrix. The intermediate steps for arriving at (46) is given in [16] and [17]. If θ^0k estimated by (45) is not on the path then it is taken as the initial point in the correction step (44). Otherwise θ^0k is considered as the next point on the path. Suppose |θ^M+1kθ^Mk|γ then we set θ^Mk=θ^k=θk+1.

8.2. Homotopy continuation algorithm

Initialize: (θ=θ0andλ=0)

do{

Increment λk+1=λk+Δλ

Update θk to θ^0k using equation (38)

if ((h(θ0k,λk+1)<threshold))

θk+1=θ0k

else

take θ0k the initial point for Newton algorithm

Update:

θ^iktoθ^k+1k using (40)

if (|θi+1kθik|threshold)

θk+1=θi+1k

else go to update: }

(Until λ=1).

Advertisement

9. Results and discussions

In simulation, two images with weak edges and two images having both weak as well as strong edges have been considered. The first original image, a liver image with ill defined edges, is shown in Fig. 3(a). In order to compute the percentage of misclassification error, the Ground Truth image, as shown in Fig. 3(b), has been constructed manually. The estimated MRF model parameters are, α= 0.005601, β = 2.34 and σ = 0.55. However σ is chosen by trial and error and is fixed at 0.5. The Percentage of Misclassification Error (PME) with respect to Ground Truth image is defined as PME = {number of misclassified pixels in all the classes}/{total number of pixels of the image}. The MAP estimates in each recursion has been obtained by our proposed hybrid algorithm [10]. The results obtained by basic MRF model is shown in Fig. 3(c), where it is observed that one of the weaker edge could be preserved while the ill defined edge adjacent to it is completely lost. In case of the CMRF model, some portions could be sharper but the adjacent ill defined edge could not be recovered as shown in Fig. 3(d). However, as seen from Fig. 3(e), the use of the proposed CCMRF model could preserve well the weak edge as well as the adjacent ill defined edges. In case of Yu 's [9] approach, the inside weak edge could not be preserved even though the outer edge could be preserved. The adjacent ill defined edge is completely lost as seen in Fig. 3(f). Thus, the proposed CCMRF model with bi-level line field with Gaussain weighted penalty function could preserve well the ill defined edges together with strong edges.

Figure 3.

(a) Liver Abscess image (468x345) (b) Ground Truth (c) MRF optimized using Hybrid (d) CMRF optimized using Hybrid (e) CCMRF optimized using Hybrid (f) Clausi’s result.

Fig.4(a) shows a cell image where the outer boundary of the cell is a strong edge while the inner portion of the cell contains weak edge or poorly defined edge. In order to compute the classification error, the corresponding ground truth image is manually constructed and is shown in Fig.4(b) and Fig.4(c) shows the result obtained with MRF model and it may be observed that the strong edges could be preserved but the weak edges could not be preserved. The poorly defined edges improved with CMRF model as shown in Fig.4(d). With CCMRF model, as observed from Fig.4(e), the outer edges of the cells could be preserved and the edges inside the different cells also have been well defined. The threshold considered for weak and strong edges are 0.91 and 0.25 respectively. The degradation process parameter is chosen to be 0.5 and the value k of the edge penalty function is chosen to be 0.2. This has also reflected in the misclassification error that is the PME is 22.72 for MRF model which reduced to 14.86 for CMRF model and further reduced to 3.11 for CCMRF model. As seen from Fig.4(f) Yu 's method preserved both weak and strong edges. The PME for Yu 's method is 6.21. It is found that the CCMRF model with bi-level line field proved to be the most effective among other methods.

Figure 4.

(a) Cell image (491x370) (b) Ground Truth (c) MRF optimized using Hybrid (d) CMRF optimized using Hybrid (e) CCMRF optimized using Hybrid (f) Clausi’s result.

In order to demonstrate the unifying modeling property of the CMRF and CCMRF model, a third example as shown in Fig. 5(a) is considered where the background has texture like attributes. The estimated MRF model parameters are α = 0.01842, β= 2.79 and σ = 0.42. As observed from Fig.5(e) that the CCMRF model could segment the image and preserved many poorly defined edges. This observation is absent in case of use of the MRF and CMRF model. Use of CCMRF model could preseve the sharp features while Yu 's method could not preseve all the weak edges. This is observed from Fig. 5(f). The percentage of misclassification error also reflect the observation. Thus, in case of all the three examples, the use of CCMRF model could segment the image and preserve both the strong as well as weak edges. This proposed model could perform better than that of Yu 's approach [9] in the context of weak edge preservation and hence misclassfication error.

Figure 5.

(a) Hand Ring (Indoor) image (303x243) (b) Ground Truth (c) MRF optimized using Hybrid (d) CMRF optimize dusing Hybrid (e) CCMRF optimize dusing Hybrid (f) Clausi’sresult.

Figure 6.

(a) MANASA SOROVER (Remote Sensing) image (500x500) (b) Ground Truth (c) MRF optimized using Hybrid (d) CMRF optimized using Hybrid (e) CCMRF optimized using Hybrid (f) Clausi’s result.

Similar observations have also been made in case of the Manasa Sorover image as shown in Fig.6(a). As observed from Fig.6(a) there are many weak edges to be preserved. In this case, the CCMRF model with bi-level linefield could preserve many weak edges together with the strong edges. This may be seen from Fig.6(e) and it can be observed from Fig.6(d) that many weak edges have been preserved even using CCMRF model. Thus, in this example the performance of CCMRF model is found to be better than that of CMRF model in the context of misclassification error.

Advertisement

Acknowledgment

Our sincere thanks for anonymous reviewers for accepting this chapter.

References

  1. 1. Y. I. Ohta, T. Kanade, T. Sakai.: “Color information for region segmentation.” Comp. Grap. Image. Process., vol 62, pp. 222-241, 1980.
  2. 2. S. Geman, D. Geman.: “Stochastic relaxation, Gibbs distributions and the Bayesian restoration of images.” IEEE. Tranaction. PAMI, vol 6, pp.721—741, 1984.
  3. 3. S.Z.Li, Markov Random Field modeling in computer vision, Springer, Berlin, 1995.
  4. 4. J. Besag: “On the statistical analysis of dirty pictures”, J.Roy. Statist.Soc.B. 62, 1986, pp.259-302
  5. 5. H. D. Cheng, X. H. Jiang, Y. Sun, Wang.J.: “Color Image Segmentation: Advances and prospects.” Pattern. Recog., vol 34, pp. 2259-2281, 2001.
  6. 6. Z. Kato, T. C. Pong, J.C>M Lee, .: “Color image segmentation and parameter estimation in a morkovian framework”. Pattern Recognition Letters, vol.22, pp309-321, 2001.
  7. 7. Y. Deng, B. S. Manjunath, : “Unsupervised Segmentation of Color-Texture Regions in Images and Video”. IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 23, no. 8, pp.800-810, 2001.
  8. 8. L.Guo, Y.M.HouandX.M.Lun, An unsupervised color image segmentation algorithm based on context information Pattern Recognition and Artificial Intelligence vol.21, pp. 82-87, 2008.
  9. 9. Qiyao Yu and David A. Clausi, “IRGS: Image Segmentation Using Edge Penalties and Region Merging, ”IEEE Transaction on Pattern Analysis and Machine Intelligence PAMI, vol. 30, no. 12, pp.2126-2139, 2008.
  10. 10. Sucheta Panda and P. K. Nanda, Constrained compound Markov Random field model with graduated penalty function for color Image Segmentation, IEEE International conference on Control, Robotics and Cybernetics(ICCRC-2011), pp.VI-126-VI-132, 2011.
  11. 11. C. Cheng, A. Koschan, C. H. Chen, D. L. Page and M. A. Abidi, Outdoor Scene Image Segmentation Based on Background Recognition and Perceptual Organization. IEEE Transactions on Image Processing, vol.21, no.3, pp.1007-1019, March2012.
  12. 12. A. K. Mishra, P. W. Fieguth, D. A. Clausi. Decoupled Active Contour (DAC) for Boundary Detection. IEEE Transactions on Pattern Analysis And Machine Intelligence, vol.33, no. 5, pp.917-930, 2011.
  13. 13. G. Scarpa, R. Gaetano, M. Haindl and J. Zerubia, Hierarchical multiple Markov chain model for unsupervised texture segmentation IEEE Transactions on Image Processing vol.18, pp.1830-1843, 2009.
  14. 14. N. Chow, J. Mallet-Paret and J. A. Yorke, Finding zeros of maps: homotopy methods that are constructive with probability one, Math.computation, vol.32, no.143, pp.887-899, 1978.
  15. 15. V. L. Stonick and S. T. Alexander, A Relationship between recursive least square update and homotopy continuation methods, IEEE Trans. Signal Processing, vol.39, no.2, pp. 530-532, 1991.
  16. 16. P. K. Nanda, K. Sunil Kumar, Sameer Gokhale and U. B. Desai, A multiresolution approach to color image restoration and parameter estimation using homotopy continuation method, Proc. IEEE Int. Conf. on Image Proc., Washington, D. C, USA, Oct. 1995.
  17. 17. P. K. Nanda, K. Sunil Kumar, Sameer Gokhale and U. B. Desai, A multiresolution approach to color image restoration and parameter estimation using homotopy continuation method, Proc. IEEE Int. Conf. on Image Processing, Washington, D.C, USA, vol 2, 45-48 Oct.1995.

Written By

P. K. Nanda and Sucheta Panda

Published: 24 October 2012