Keywords

1 Introduction

Texture analysis is an active area of research in image processing. The main reason for the importance of texture analysis is that its applications are pervasive in many fields of study, such as biometric identification [1], remote sensing [2], science and technology [3], biology [4], medicine [5], and visual arts [6]. However, the concept of texture in images is inherently imprecise, making it difficult and limited for the mathematical formulation of its complex property. Furthermore, noise is well known as a major factor that hinders the performance of extracting effective texture features for image classification, motivating many efforts to develop techniques for removing or robustly working with noise in images [7,8,9,10].

In general, noise in signals is undesirable for image analysis, but there are circumstances when some certain amount of noise can be useful, for example, to prevent discretization artifacts in color banding or posterization of an image [11]. The preserve of film grain noise can also help enhance the subjective perception of sharpness in images, known as acutance in photography, although it degrades the signal-to-noise ratio [12]. The intentional inclusion of noise in processing digital audio, image, and video data is called dither [13, 14].

Another useful aspect of noise for texture analysis is reported in this paper. The motivation is mainly based on the concept of geostatistics that employs random models to characterize spatial attributes of natural phenomena [15]. In statistics, numerical methods are developed to deal with variables or random variables. However, many variables of natural phenomena that vary in space and/or time may not be all random. Some variables may be totally deterministic, and some may take place somewhere between randomness and. Thus, the concept of a regionalized variable as a variable distributed in space is introduced to capture a behavior that is characterized with both of randomness and determinism [15]. In this regard, a regionalized variable is deterministic within a spatial domain, but beyond that it behaves as a random variable. The semi-variogram of geostatistics is formulated to quantify the behavior of a regionalized variable [16]. In this sense, adding noise at some certain levels into a texture space can give rise to discriminative power in the quantification of the random part of a regionalized variable of the image intensity by means of the semi-variogram. Due to the limitation of the semi-variogram for expressing the spatial continuity of a regionalized variable with experimental data, a spectral distortion measure is adopted to measure the dissimilarity of regionalized variables for texture retrieval.

2 Regionalized Variables in an Image

Since the intensities of pixels are variables distributed in space, they can be modeled as regionalized variables, each of which is considered as a single realization of a random function. In one sense, this regionalized variable is not related to its neighboring variables (pixels). In other sense, this variable has a spatial structure, depending on the distance separating the pixels. Thus, with the combination of the random and structured properties of a regionalized variable in a single function, the spatial variability of an image can be described on the basis of the spatial structure of these variables [17].

Without loss of generality, let Z(x) be a regionalized variable, which is a realization of a random function Z, x and h be a spatial location and a lag distance in the sampling space, respectively. The variogram of the random function is defined as [16]

$$\begin{aligned} 2 \gamma (h) = Var[Z(x)-Z(x+h)], \end{aligned}$$
(1)

where \(\gamma (h)\) is the semi-variogram of the random function. This definition of the variogram, \(2\gamma (h)\), or semi-variogram, \(\gamma (h)\), assumes that the random function changes within the space, but \(\gamma (h)\) is independent of spatial location and depends only on the distance of the pair of the considered variates. To simplify technical jargon, the semi-variogram is now referred to as the variogram, unless mathematical expression requires a precise definition.

Based on Eq. (1), the variogram is equivalent to

$$\begin{aligned} \gamma (h) = \frac{1}{2} E\left[ \{Z(x) - Z(x+h)\}^2\right] \end{aligned}$$
(2)

Let \(Z(x_i)\), \(i=1, 2, \dots , n\), be a sampling of size n, the unbiased estimator for the variogram, which is called the experimental variogram, of the random function is expressed as

$$\begin{aligned} \gamma (h) = \frac{1}{2N(h)} \sum _{i=1}^{N(h)} \left[ Z(x_i) - Z(x_i+h) \right] ^2, \end{aligned}$$
(3)

where N(h) is the number of pairs of variables separated by distance h.

The formulation of Eq. (3) is based on the assumption that the spatial autocorrelation structure is isotropic. This means that the semi-variogram depends only on the magnitude of the lag (h). When the spatial autocorrelation pattern is not the same in different directions in the sampling space, an anisotropic semi-variogram should be used to accommodate these differences. There are two types of anisotropy of the semi-variogram: geometric and zonal [18]. The geometric anisotropy occurs when the range of the semi-variogram varies with different directions. The zonal anisotropy takes place when both the range and sill of the semi-variogram change with different directions. The sill, which is the upper bound of the semi-variogram, represents the variance of the random field, whereas the range, at which the semi-variogram reaches the sill, indicates the distance at which data are no longer autocorrelated. These two parameters of the semi-variogram expressed by the spherical (theoretical) function [16], \(\gamma ^T(h)\), are shown in Fig. 1.

Fig. 1.
figure 1

Theoretical semi-variogram using spherical model with \(s=1\) (sill) and \(g=20\) (range)

Full size image

With two coordinates for a 2D case, where \({\mathbf h} = (h_1, h_2)\), then a model for the geometric anisotropic semi-variogram can be defined as [18]

$$\begin{aligned} \gamma (r) = \gamma \left( \sqrt{{\mathbf h}^T {\mathbf B} {\mathbf h}} \right) , \end{aligned}$$
(4)

where \({\mathbf B} = {\mathbf Q}^T {\mathbf \Lambda } {\mathbf Q}\), and \({\mathbf Q}\) is the transformation matrix:

$$\begin{aligned} {\mathbf Q} = \left[ \begin{array}{cc} \cos \theta &{} \sin \theta \\ -\sin \theta &{} \cos \theta \end{array} \right] , \end{aligned}$$
(5)

where \(\theta \) is the rotation angle, and \({\mathbf \Lambda }\) is the diagonal matrix of eigenvalues:

$$\begin{aligned} {\mathbf \Lambda } = \left[ \begin{array}{cc} \lambda _1 &{} 0 \\ 0 &{} \lambda _2 \end{array} \right] . \end{aligned}$$
(6)

A model for the zonal anisotropic semi-variogram is defined as [18]

$$\begin{aligned} \gamma ({\mathbf h}) = \gamma _1({\mathbf h}) + \gamma _2({\mathbf h}), \end{aligned}$$
(7)

where \(\gamma _1({\mathbf h})\) is the isotropic semi-variogram in one direction whose sill is much larger than the sill produced in the other direction, and \(\gamma _2({\mathbf h})\) is the geometric anisotropic semi-variogram.

3 Measuring Dissimilarity of Regionalized Variables

Given an experimental variance at lag h, \(\gamma (h)\) has recently been proposed to be approximated as a linear combination of the past p variances [19]

$$\begin{aligned} \tilde{\gamma }(h) = - \sum _{i=1}^p a_i \gamma (h-i) \end{aligned}$$
(8)

where \(a_i\), \(i=1, \dots , p\) are the linear predictive coding (LPC) coefficients [20], and to be optimally determined as follows.

The error between \(\tilde{\gamma }(h)\) and \(\gamma (h)\) is given by

$$\begin{aligned} e(h) = \gamma (h) + \sum _{i=1}^p a_i \gamma (h-i) \end{aligned}$$
(9)

By minimizing the sum of squared errors, the pole parameters \(\{a_i\}\) of the LPC model can be determined as follows.

$$\begin{aligned} {\mathbf a} = - {\mathbf R}^{-1} \, {\mathbf r} \end{aligned}$$
(10)

where \({\mathbf a}\) is a \(p \times 1\) vector of the LPC coefficients, \({\mathbf R}\) is a \(p \times p\) autocorrelation matrix, and \({\mathbf r}\) is a \(p \times 1\) autocorrelation vector whose elements are defined as

$$\begin{aligned} r_i = \sum _{h=0}^N \gamma (h) \gamma (h+i), \, i = 1, \dots , p. \end{aligned}$$
(11)

Let \(S(\omega )\) and \(S'(\omega )\) be the spectral density functions of the semi-variograms \(\gamma (h)\) and \(\gamma '(h)\), respectively, where \(\omega \) is a normalized frequency ranging from -\(\pi \) to \(\pi \). The spectral density \(S(\omega )\) is defined as [20]

$$\begin{aligned} S(\omega ) = \frac{\sigma ^2}{|A|^2}, \end{aligned}$$
(12)

where \(\sigma ^2 = {\mathbf a}^T{\mathbf R} {\mathbf a}\), and \(A = 1 + a_1 e^{-i\omega } + \dots + a_p e^{-ip\omega }\).

The log-likelihood-ratio (LLR) distortion measure between \(S(\omega )\) and \(S'(\omega )\), denoted as \(D_{LLR}(S,S')\), is defined as [21]

$$\begin{aligned} D_{LLR}(S,S') = \log \frac{{\mathbf a}'^T {\mathbf R} {\mathbf a}'}{{\mathbf a}^T {\mathbf R} {\mathbf a}}, \end{aligned}$$
(13)

where \({\mathbf a}'\) is the vector of the LPC coefficients of \(S'\).

4 Experiments

The proposed method was tested using four subsets of the Brodatz database [22] to represent four manually-labeled types of texture [23]: (1) fine-periodic, (2) fine-aperiodic, (3) coarse-periodic, and (4) coarse-aperiodic. Each subset consists of 90 images of \(215 \times 215\) pixels, which are produced by dividing each of the 9 corresponding original Brodatz images into 9 non-overlapping smaller samples. The Brodatz indices of the fine-periodic texture are: D3, D6, D14, D17, D21, D34, D36, D38, D49, and D52. The image indices of the fine-aperiodic texture are: D4, D9, D16, D19, D24, D26, D28, D29, D32, and D39. For the coarse-periodic texture, the ten images are: D1, D8, D10, D11, D18, D20, D22, D25, D35, and D47. For the coarse-aperiodic texture, the ten images are: D2, D5, D7, D12, D13, D15, D23, D27, D30, and D31.

The isotropic model of the semi-variogram was used to compute the semi-variances for 30 lags of all the images, because the ranges of the images are almost the same, the use of the isotropic model has been found suitable for extracting image features, and the computation of the anisotropic variogram is very time-consuming, thus, it is not possible for large images.

The images were degraded with white Gaussian noise (additive noise), speckle noise (multiplicative noise), and salt & pepper noise (impulsive noise). Figure 2 shows a sub-image of D30 representing a coarse-aperiodic texture degraded with white Gaussian noise of zero mean and 0.1 variance, speckle noise with zero mean and 0.1 variance, and salt & pepper noise with 0.1 noise density.

Fig. 2.
figure 2

Brodatz coarse-aperiodic D30 (a), and corresponding degraded images with white Gaussian noise with zero mean and 0.1 variance (b), speckle noise with zero mean and 0.1 variance (c), and salt & pepper noise with 0.1 noise density (d).

Full size image
Fig. 3.
figure 3

Retrieval rates of four Brodatz image subsets degraded with various noise models for four texture categories: fine-periodic (a), fine-aperiodic (b), coarse-periodic (c), and coarse-aperiodic (d); where the zero noise level indicates the corresponding retrieval rate of undegraded images.

Full size image
Fig. 4.
figure 4

Retrieval rates of four Brodatz subsets degraded with Gaussian noise with 0.01 variance and variable means, where the zero Gaussian noise mean indicates the corresponding retrieval rate of undegraded images.

Full size image

For the matching of pattern similarity between two images, \(D_{LLR}\) was used with \(p=20\), to compute the spectral distortion of the semi-variograms of the two corresponding images. The retrieval was carried out for each of the 360 (\(10 \times 9 \times 4\)) degraded images by searching for 8 images of the same texture among the most k similar images, where \(k=21\) in this case. Figure 3 shows the plots of the retrieval rates of the four Brodatz image subsets degraded with three noise models and a variety of noise levels.

For the fine-periodic texture (Fig. 3(a)), the retrieval rates of images degraded with salt & pepper noise maintain the same rate in comparison with the undegraded images up to the noise level of 0.06, whereas 0.2 and 0.3 for Gaussian noise and speckle noise, respectively. For the fine-aperiodic texture, Fig. 3(b) shows the consistently high performance of the images added with speckle noise, whereas the retrieval rate of the images degraded with Gaussian noise with 0.01 variance is higher than that of the undegraded images, and the rate using the inclusion of salt & pepper noise at 0.01 density level is slightly lower than the retrieving performance without noise. For the coarse-periodic texture shown in Fig. 3(c), the retrieval rates are higher for images degraded at the noise level of 0.01 using all the three noise models than the retrieval rate for undegraded images, where speckle noise gives the highest result. The same results also apply to the coarse-aperiodic (Fig. 3(d)), where Gaussian noise has the best performance. Keeping the Gaussian noise variance constant at 0.01, while varying the Gaussian noise mean from 0.01 to 0.1, with the exception of the fine-aperiodic texture, Fig. 4 shows the high performance of the retrieval task by adding Gaussian noise to the images of the other three types of texture.

5 Conclusion

Texture analysis using the spectral distortion measure of semi-variograms can be enhanced by adding noise at certain levels to the images. The set of 40 Brodatz images were selected in this study because they fairly represent the 4 types of texture. The addition of different types of noise to other types of texture is worth pursuing further investigation. Furthermore, in this paper, only empirical evidence was provided to support the argument of the useful application of noise to texture, theoretical verification of the power of geostatistics in finding structures in random-biased texture is therefore important to follow.