Neighbouring Color Dependence Matrix for Image Analysis : Application to homogeneous and heterogeneous areas detection and characterization

Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 Neighbouring Color Dependence Matrix for Image Analysis : Application to homogeneous and heterogeneous areas detection and characterization B. Jacquin∗and A. Smolarz∗ ∗University of Technology of Troyes, ICD Laboratory FRE CNRS 2848, BP 2060 F10010, TROYES CEDEX, FRANCE Received 12 may 2006; revised 4 september 2006; accepted 26 February 2008 Abstract A new method for color texture characterization and color texture region detection is presented. This method, which we will name NCDM (Neighbouring Color Dependence Matrices), is the extension to color textures of the NGLDM (Neighbouring Gray Level Dependence Matrices) introduced by Sun et al. [1] and completed by Berry et al. [2]. This approach consists in estimating the dependences of colors between a pixel and its neighbours. We propose two steps: a color areas classification in two classes followed by the characterization of the detected areas. In the first step, we compute the NCDM with an isotropic neighbourhood. The structure of the isotropic NCD distribution allow us to separate the pixels of a color composite image into two classes, which correspond respectively to homogeneous and heterogeneous regions in the image. We then consider that the heterogeneous regions are potentially textured regions and in the second step we propose to compute the NCDM with anisotropic neighbourhoods corresponding to the eight principal directions. To seek the dominant directions in a color texture, a measure of spatial dependence between a pixel and its neighbours is computed by way of a chi-square test. This measure is based on the fit of the NGLD and NCD distribution with a binomial model under independence hypothesis. The variations of the colors are computed in uniform perceptual color spaces. We have chosen the color space ”L1 norm” introduced by Angulo and Serra [3]. Key Words: Color space, Color Texture, Anisotropy, NGLDM, Color Image Segmentation, Chi-square Test. 1 Introduction We propose a new method for the analysis of color images which is of important use in content-based image retrieval (CBIR) systems, composed of two principal steps. The first step consists in separating a color image into two classes. In order to do that, we classify all the pixels that satisfy a criterion of homogeneity (in terms of color) as belonging to a ”homogeneous class”. Let us specify that class relates to an area because our criterion of homogeneity is a spatial criterion with a local measure in the neighbourhood of a pixel. The non-labeled pixels do not satisfy the criterion of homogeneity that we propose, and these pixels can belong to textured areas, to noisy areas, or to the boundaries between the areas. The second step consists in extracting a set of features on these pixels in order to detect and to characterize the heterogeneous areas. The features are based on a measure Correspondence to: <[email protected]> Recommended for acceptance by Jean-Marc Ogier ELCVIA ISSN:1577-5097 Published by Computer Vision Center / Universitat Aut` onoma de Barcelona, Barcelona, Spain 2B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 of anisotropy in the eight principal directions. Our method is an extension of the Neighbouring Gray Level Dependence Matrices (NGLDM, Sun & Wee in 1983 [1] and Berry & Goutsias in 1991 [2]) to color images, which we will call : Neighbouring Color Dependence Matrices (NCDM). These matrices evaluate the degree and the nature of the local dependence in the neighbohood of a pixel. The definition of a NCDM matrix will be detailed in section 3. To evaluate the efficiency of the proposed method, we will use a publicly image database and we will give the number of correct classified pixels. Initially, we will review the main works that deal with the definition of the color texture (section2). In section 3, we recall the definition and the notations used with NGLD matrices. We also present the measures of homogeneity and heterogeneity that we propose for the first step of classification. Then we define the measure of anisotropy and the following set of features used for the second step of detection and characterization of texture. Finally, this section is completed with the extension of the definition of NGLD matrices to NCD matrices and we discuss the consequences of applying the methodology to color images. In section 4, we present the first step of separation of the pixels in two classes corresponding to the homogeneous and heterogeneous areas. Some examples will illustrate this first step which is computed in RGB color space. The second aspect of our paper is presented in section 5. By considering that the heterogeneous areas are potentially textured areas, we calculate the NCD matrices associated with anisotropic neighbourhoods corresponding to the eight principal directions. Some examples will illustrate this second step which is computed in ”L1 norm” color space ([3] and [7]). Finally, we will present our conclusion and the prospects for future work. 2 Previous works We describe here works on color textures. Firstly, we will discuss the perceptual properties that can be observed when a texture is examined. Rao and Lohse [9] reported the results of an experiment in which the aim was to find the high level characteristics of the texture for human vision. Twenty people were asked to independently classify 30 textures of the Brodatz album [12]. The authors identified the most significant characteristics as the regularity (a pattern which is repeated), followed by the orientation and the complexity. Mojsilovic et al. [13] presented the results of a similar experiment with color textures. They asked 28 subjects to specify numerically, on a scale from 0 to 100, the similarities between each possible combination of 20 color textures extracted from a set of fabrics. After a statistical analysis, they identified five significant perceptual characteristics. Two characteristics were specific to the color: the dominant color and the purity of the color. The three other characteristics corresponded to the characteristics found by Rao and Lohse: directionality and orientation, regularity, complexity and density. The NCD matrices, which are presented in the next section, implicitly provide these features : the dominant color, the orientation and the structure of the textures. Angulo, Hanbury and Serra in [5], [4] and [3] have worked with oriented textures for studies on the automatic classification of wood according to two criteria: the color and the texture. They sought the dominant directions in the neighbours of each point of the image. Serra and Hanbury used the algorithm proposed by Rao ([10] and [11] ) to find these directions. The first step of Rao’s algorithm ( [10] and [11] ) applies a Gaussian filter to the texture. It is necessary to calculate an angle for each pixel with the images for horizontal and vertical gradient. The dominant angle is calculated in the neighbourhoods in order to produce an angle image. Each pixel in the angle image corresponds to the dominant direction in a neighbourhood. Germain et al. [14] worked on the characterization of the anisotropy of the textured images. In their works, they reported that the dominant direction or the anisotropy in the textures strongly depends on the observation scale. They introduced a new operator for the estimation of the dominant direction, the Directional Mean Vector (DMV ), which works only on gray level textures. With the NCD matrices, different sizes of the neighbourhood enable a multi-scale analysis to be used in order to characterize the anisotropy in the color textures. Shu and Jain [15] presented a method based on the properties of vector fields for the estimation of a set of symbolic descriptors from linear orientation fields. The authors presented results of experiments on real oriented B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 3 textures. Almeida [16] presented a Markov random field (MRF) model for digital images that were able to represent anisotropic textures with arbitrary orientations. Works relating to the characterization of textures in gray levels with NGLD matrices were proposed by Smolarz in 2003 in [6]. An extension to color textures of these works is presented in section 5 of the current work. Fontaine et al. [17] presented tools for spatial and colorimetric analysis that take into account the global and local colorimetric information of an image, such as the color co-occurrence matrices and the color correlograms initially introduced by Huang et al. [18]. 3 Presentation of the NGLDM and NCDM In order to clearly present the criteria that we have defined on the NCDM, we will recall some notations and definitions related to the NGLDM. 3.1 Notations used for the NGLDM NGMaximum number of gray levels. G(`, c)Random variable that represent the gray level of the pixel (`, c).G∈ {0,1, . . . NG−1} dParameter equal to the maximum distance (in pixels) between the central pixel and its neighbours. VdThe set of neighbours, called the neighbourhood. NdThe number of neighbours. Nd= #{Vd} aThe threshold of the difference between the gray levels of the central pixel and one of its neighbours. SRandom variable that represents the number of neighbours which have a gray level ”close to” (as defined by parameter a) that of the central pixel. S∈ {0,1, . . . Nd} In the original version of Sun et al. [1], the parameters dand aof the NGLDM are fixed by the user. In the current application, the parameter awill be fixed automatically for the NCD matrices (see section 4), which reduces the need for user input, as well as providing greater objectivity. 3.2 Definition of the NGLDM In order to simplify the notations, the NGLDM will be noted as Q. The general term of Qis given by the expression: Qd,a(g, s) = # {(`, c)|G(`, c) = gand # [(i, j)∈Vd| |G(`, c)−G(i, j)| ≤ a] = s}(1) Hence the NGLDM is a matrix of size NG×(Nd+ 1) and Qd,a(g, s)is the number of pixels in an image with gray level gthat have exactly sneighbours with gray level in the interval [g−a , g +a]. The matrix Qallows us to calculate the estimated joint distribution of the random variables Gand S, as shown below: ˆ Pd,a [G=g , S =s] = Qd,a(g, s) PNG−1 u=0 PNd v=0 Qd,a(u, v)(2) 3.3 Homogeneity, heterogeneity and NGLDM Firstly, in order to identify the homogeneity, we make the same assumptions as Deng et al. [19]: each image contains a set of approximately homogeneous color-texture regions. The color information in each image region can be represented by a set of a few quantized colors. The colors between two neighbouring regions are distinguishable. A region is called heterogeneous region if each pixel and its surrounding neighbours have 4B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 different colors. If we consider only the first columns of Q, it is possible to find the gray levels that best characterize the heterogeneous regions. Conversely, if we consider the last columns of Q, it is possible to find the gray levels that best characterize the homogeneous regions. For example, consider the gray level g: if for n1Ndwe have a high value of Pn1 v=0 ˆ Pd,a [G=g , S =v] we can conclude that the gray level glies in a heterogeneous region. If for 0n2< Ndwe have a high value of PNd v=n2 ˆ Pd,a [G=g , S =v]we can conclude that the gray level glies in a homogeneous region. Our approach to distinguish between homogeneous and textured regions is based on these properties of the NGLDM. The method will be detailed and illustrated in section 4. 3.4 The binomial model and the anisotropy measure According to the previous formula (2), we can define the following conditional probability estimation: ˆ Pd,a [S=s|G=g] = Qd,a(g, s) PNd v=0 Qd,a(g, v)(3) Under the non realistic (especially in the case of texture modeling) hypothesis that the random variables {G(`, c)|(`, c)∈Vd}are independant, the distribution ˆ Pd,a [S=s|G=g]estimated by (3) is a binomial distribution with parameters Ndand p=P[g−a≤G≤g+a]. We have choosen to use the well known Chi-Square distance to test the independence hypothesis. The resulting criterion is computed as follows: 1. For each gray level g, compute the Chi-Square distance : E(g) = Nd X s=0 Qd,a(g, s)−Qd,a(g, ·)ˆ Pd,a [s|g]2 Qd,a(g, ·)ˆ Pd,a [s|g]where Qd,a(g, ·) = Nd X s=0 Qd,a(g, s) 2. For a fixed type one error α(typically 5%), calculate the difference D(g) = E(g)−χ2 αbetween the ChiSquare distance and the test threshold at level α. 3. Finally compute the mean value D=1 NG NG−1 X g=0 D(g)(4) If Dhas a negative value, we can accept the independence hypothesis. In the case of a positive value for D, the greater the value of D, the stronger the dependence hypothesis. To take the anisotropy of textures into account, we can compute many NGLDM defined on anisotropic neighbourhoods, which are simply lines. Consider, for example, the four principal directions 0◦,45◦,90◦and 135◦. We obtain four values D0, D45, D90 and D135 of the criterion defined by (4). {D0, D45, D90, D135}constitutes a set of features that can be used to measure the anisotropy of a texture. In order to have normalized features, we propose to use the set {f0, f45, f90, f135}defined by: fθi=Dθi Pj Dθj  and so fθi∈[−1,1] ∀θi= 0◦,45◦,90◦,135◦(5) In the case of independence (for example a uniform or gaussian noise), all the features fθiwill have the theoretical value −0.25, whereas in the case of an isotropic texture, all the features will have values close to 0.25. These properties have been tested, with full details available in [6] B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 5 3.5 Definition of the NCDM We will now investigate how to define a Neighbouring Dependence Matrix based on color images. If we refer to formula (1), the main problem comes from the expression |G(`, c)−G(i, j)| ≤ a. The question is how can the difference between two colors be calculated ? The answer that we propose is, according to a pre-defined color space, to use a distance between two colors that is called Perceptual Difference P D (see [4]). The Perceptual Difference can be defined as follows. We consider two color spaces, the cartesian space RGB and the cylindrical coordinate space HBS. The gray level random variable G(`, c)in formula (1) will be substituted by a random color variable C(`, c) = [R(`, c), G(`, c), B(`, c)] in the cartesian color space and C(`, c)=[H(`, c), B(`, c), S(`, c)] in the cylindrical coordinate space. First of all, let us introduce some simple notations that will be used afterwards. Let C1=C(`1, c1)and C2=C(`2, c2)be the colors of two pixels (`1, c1)and (`2, c2)then, P D(C1, C2)is defined by: P D(C1, C2) = p(R1−R2)2+ (G1−G2)2+ (B1−B2)2(RGB space)(6) P D(C1, C2) = q(B1−B2)2+S2 1+S2 2−2S1S2cos(H1−H2) (HBS space)(7) According to the definitions given by the formulas (6) and (7), we can now define a new matrix Qwhich we call NCDM with the general term: Qd,a(h, s) = # {(`, c)|C(`, c) = hand # [(i, j)∈Vd|P D(C(`, c), C(i, j)) ≤a] = s}(8) It should be noted that hin formula (8), represents a color triplet rather than a scalar as is the case for gin formula (1). Let us clarify our model for the measure of the anisotropy in a color image: ˆ Pd,a [S=s|C=h]is the estimated distribution, conditionally with the color of a pixel, of the number sof its neighbours having a close color (the perceptual difference P D is lower than a). Always under the assumption that the pixels belonging in the neighbourhood are independent, the probability ˆ Pd,a [S=s|C=h]can be modeled by the binomial distribution with parameters Ndand p=P[DP ≤a]. The test is carried out by computing the Chi-2 distance for each color hwith pdependent on h[7]. Then, we compute the difference D(h)between the Chi-2 distance and the decision threshold χ2 α, obtained for a fixed type one error, in this case α=5%. In order to take the color distribution into account, we calculate the mean of the differences D(h)balanced by the frequency (see formula (9)) instead of the arithmetic mean suggested with the formula (5). We compute this mean for each direction θi. Dθi= Nh−1 X h=0 D(h)f(h)where f(h) = Nd X v=0 ˆ Pd,a [C=h, S =v](9) where f(h)represents the frequency of the color h. Finally, as for gray level images, we compute the standardized features according to formula (5). Negative values of fθiindicate that there is no dependence between the pixels in the selected direction. The larger the positive coefficients fθi, the greater the dependence between the pixels. 3.6 Color reduction and consequences on NCD matrices The size of a NCD matrix is Nh×(Nd+1) where Nhis the number of colors in the image. A composite image can have several thousands of colors, and as for the co-occurrence matrices, we need to reduce the number 6B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 of colors in order to avoid the estimation problems of NCD distribution. This quantification step of colors is critical for the separation step of the pixels into two classes corresponding to the homogeneous and the heterogeneous areas. Indeed, the variations of contrast and hue must not be amplified during this phase, particulary for homogeneous areas and spectrally homogeneous regions (reflections, shadows). For optimal results, the quantification step proceeds after having changed color space in order to exploit the colorimetric properties of the new color space. After several tests of color reduction methods [20], we have retained the minimum variance quantization algorithm, which creates a new pallete of colors. After this algorithm, there is a spatial quantification step to smooth the contrast areas. With this algorithm, we can choose the number of colors in the new palette. In the separation step of the pixels into homogeneous areas and heterogeneous areas, we reduce the composite image to 16 colors in order to have a good compromise between the estimation constraints of NCDM and the degradation of the homogeneous areas. This choice implies the assumption that there is, at most, 16 different homogeneous areas in a composite image, which seems to be a reasonable assumption. Let us specify that the number of colors is a parameter. We also tested the two steps by reducing the images to 32 colors. The results are similar and are presented on figure 3. For the next step of anisotropy characterization in the potentially textured areas, we reduce the remaining areas to 64 colors, because the texture characterization requires a higher level of details in the colors. In order to work on estimated probabilities, NCD matrices are normalized according to formula (2) and with the notations of the formula (8). The construction of an NCD matrix depends on two parameters (see section 3): the parameter d, which gives the size of the neighbourhood and the parameter a, which is the threshold of the difference in color between two neighbouring pixels. The avalue changes according to the analyzed image and thus, we propose to determine it automatically for each image. After color reduction of an image, we calculate the histogram and then, the avalue is fixed automatically by seeking the mean of the perceptual distances (P D) for all the color couples. This mean is computed as follows. a= Nh−1 P u=0 Nh−1 P v=0,v6=u P D(u, v)f(u) Nh(Nh−1) (10) In the two steps of this paper, we use a very interesting property of NCD matrices : computing a NCD matrix on an image which contains thousands of colors with an adequate value for ais equivalent to computing a NCD matrix with a reduced color number and a larger value for a. Indeed, increasing the avalue is equivalent to summing the information of several lines of NCD matrix, and therefore, merging several colors. Thus, NCD matrices are compatible with a color reduction step. 4 Separation step between homogeneous and heterogeneous regions As indicated previously, for this step, we reduce the analyzed image to 16 colors. We choose an isotropic neighbourhood with a size d= 2 (Nd= 24 neighbours). Thus, the NCD matrix has a size of 16 ×25. Let us take for example the ”Mandrill” image (presented on figure 2.a) reduced to 16 colors. This image can be seen as the fusion, or in other words the union, of 16 color labeled layers I(h). Each color labeled layer I(h) corresponds to 1 color among the 16 colors (see figure 1). These color labeled layers I(h)are complementary, and form a partition on the set of the pixels of the initial image. Our goal is to sort these color labeled layers I(h)from the most homogeneous to the most heterogeneous. In section 3.3, we introduced two thresholds n1and n2. We chose here n1= 4 and n2= 20, which means that the areas where a pixel has at most four neighbours of close color (≈20% of Nd) are labeled as heterogeneous and the areas where a pixel has at least twenty neighbours of close color (≈80% of Nd) are labeled homogeneous. Thus, for a given color h, if ˆ Pd,a [C=h , S =s]has a maximum in the left hand part of the matrix (s≤n1), then his mainly present in the heterogeneous areas of the image. Conversely, if B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 7 ˆ Pd,a [C=h , S =s]has a maximum in the right hand part of the matrix (s≥n2), then his mainly present in the homogeneous areas of the image. r(h)=0.53 r(h)=0.31 r(h)=0.19 r(h)=0.15 r(h)=0.13 r(h)=-0.05 r(h)=-0.10 r(h)=-0.15 r(h)=-0.15 r(h)=-0.20 r(h)=-0.29 r(h)=-0.34 r(h)=-0.36 r(h)=-0.99 r(h)=-0.40 r(h)=-0.45 Figure 1: Layers of a ”Mandrill” sorted from the most homogeneous one to the most heterogeneous one. The r(h)values for each color labeled layer is displayed above it Then, for each color h, that is to say for each color labeled layer I(h), we compute the estimated probability het(h)that a color belongs to a heterogeneous area and the estimated probability hom(h)that a color belongs to a homogeneous area according to the formulas (11) and (12). het(h) = ˆ Pd,a [C=h, S ≤n1] = n1 X v=0 ˆ Pd,a [C=h, S =v](11) hom(h) = ˆ Pd,a [C=h, S ≥n2] = Nd X v=n2 ˆ Pd,a [C=h, S =v](12) In order to work with a relative and standard measure, we compute a ratio r(h)as follows. −1≤r(h) = hom(h)−het(h) f(h)≤+1 where f(h)represents the frequency of the color h(13) The color labeled layers are then sorted according to the values r(h)with decreasing homogeneity (figure 1). It should be noted that the size of the homogeneous or heterogeneous areas (that is to say the number of pixels for each color labeled layer) does not influence this classification because the terms in the NCD matrix are normalized (see formula (2)). (a) original image (b) homogeneous areas (c) heterogeneous areas (d) profile of r(h) mean of r(h)=-0.11 Figure 2: An image of a ”Mandrill” with 16 colors. The homogeneous and heterogeneous areas in (b) and (c) result from the threshold aobtained from figure 1 8B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 (a) homogeneous areas (b) heterogeneous areas (c) profile of r(h) mean of r(h)=-0.18 Figure 3: An image of a ”Mandrill” with 32 colors. The homogeneous and heterogeneous areas in (a) and (b) were obtained with a threshold aequal to −0.18 From now on, we have to find a threshold Tfor r(h)such that the fusion of the color labeled layers I(h) with r(h)≤Tforms the set of the potentially textured areas and the fusion of the color labeled layers I(h) with r(h)> T forms the set of the homogeneous areas of the image. Experimentally, we have noticed that the arithmetic mean value of r(h)seams to be a good threshold to separate between homogeneous and heterogeneous areas. A profile of r(h)appears on figure 2.d (his on the X-axis and r(h)on the Y-axis). We can see that the mean value is located towards the zero-crossing of this curve. So an other possible choice of the threshold Tcan be the value hof zero-crossing. The results of our experiments indicated that the mean value of r(h)gives better classification results between hom(h)and het(h)on the test set. Remembering that the color labeled layers have been sorted in decreasing order (according to r(h)), the threshold Tcorresponds to the color hof the less homogeneous layer I(h), i.e. the layer with the value of r(h) closest to T. For the example of figure 1, the arithmetic mean value of r(h)is −0.11. Therefore, the color labeled layer corresponding to the threshold will be the one with r(h) = −0.10. This selected layer will be merged with all the color labeled layers more homogeneous than it (indicated by an orange stripe above them in figure 1), to form the set of homogeneous areas of the image ”Mandrill” (figure 2.b). The other color labeled layers are merged, in order to form the set of heterogeneous (potentially textured) areas, which correspond to the fur of the Mandrill (figure 2.c). The final result is a partition of the original image into two classes : the homogeneous and the heterogeneous areas. To illustrate the separation step of the pixels into two classes, we present, in figure 4, the results obtained on an image of a rabbit (the fur is well isolated), a carpark (the trees are well isolated), a golf course (the vegetation and its reflection in the water are identified as heterogeneous), some flowers (the flowers are well isolated) and a human face (the skin is well isolated). We can see that an heterogeneous area can be composed of a frontier between two homogeneous areas or a degraded gray-level. Let us note that this step of separation is very precise since the computation is made at a pixel level, thus enabling very fine edges to be obtained (see, for example, the tomato in the image of the rabbit). This first step, from the loading of the image, thought the analysis until the creation of the two classes, requires less than one second of computing time for images of 512 ×512 pixels. In order to evaluate in a quantitative way the efficiency of our approach, we have generated four synthetic images with known homogeneous and heterogeneous regions (see figure 5). The heterogeneous areas are composed of textures downloaded from the image database available at http://textures.forrest.cz/. In the four examples the background forms the homogeneous areas. The performance criterion is the rate of correct classified pixels calculated after the detection step of the two classes. For the first three images (figures 5.a, 5.d and 5.g), the percentages of correct classified pixels are respectively 94.4%,99.7% and 97.7%. These values, B. Jacquin and A. Smolarz / Electronic Letters on Computer Vision and Image Analysis 7(1):1-15, 2008 9 original image homogeneous areas heterogeneous areas profile of r(h) original image homogeneous areas heterogeneous areas profile of r(h) original image homogeneous areas heterogeneous areas profile of r(h) original image homogeneous areas heterogeneous areas profile of r(h) original image homogeneous areas heterogeneous areas profile of r(h) mean of r(h)=0.02 mean of r(h)=0.07 mean of r(h)=-0.07 mean of r(h)=0.08 mean of r(h)=0.03 Figure 4: Separation of the pixels into two classes : homogeneous areas and heterogeneous areas, and the profiles of r(h)