Patterns Simulations Using Gibbs/MRF Auto-Poisson Models

Abstract

Pattern analysis is the process where characteristics of big data can be recognized using specific methods. Recognition of the data, especially images, can be achieved by applying spatial models, explaining the neighborhood structure of the patterns. These models can be introduced by Markov random field (MRF) models where conditional distribution of the pixels may be defined by a specific distribution. Various spatial models could be introduced, explaining the real patterns of the data; one class of these models is based on the Poisson distribution, called auto-Poisson models. The main advantage of these models is the consideration of the local characteristics of the image. Based on the local analysis, various patterns can be introduced and models that better explain the real data can be estimated, using advanced statistical techniques like Monte Carlo Markov Chains methods. These methods are based on simulations where the proposed distribution must converge to the original (final) one. In this work, an analysis of a MRF model under Poisson distribution would be defined and simulations would be illustrated based on Monte Carlo Markov Chains (MCMC) process like Gibbs sampler. Results would be illustrated using simulated and real patterns data.

1. Introduction

Pattern analysis is a class of methods that are used to recognized regular patterns in big data, like images. Considering these methods, models must be introduced, and analysis of the modeling process must be performed. Especially for the images, the spatial structure of the pixels is an important measure to explain the spatiality of the image.

Spatial structure of images can be analyzed based on specific models, where neighborhood structures are taking into consideration. Interactions between regions at different scales are characterized by their local dynamics, and the emergent spatial patterns are the outcome of different processes. Local dynamics could be explained by the local characteristics of the image with the main result being the construction of the images under homogeneous regions. The homogeneity of the image could be defined by the spatial pattern of the image, which is explained with way better presentation of real data. Under these local characteristics, special models considering the neighborhood structure of the image could established specific conditional models defined as Markov Random fields (MRF) models. These models are defined by Markov Random fields (MRF) models, explaining the spatial structure of the images by a conditional distribution of the pixels, defined as auto-models. These models under appropriate structure patterns have the ability to explain real patterns. Under the conditional distribution, various models could be introduced. In our work, the specific class of data that can be analyzed as best as it can is based on the Poisson distribution, defined as auto-Poisson models.

Referring to conditional probabilities, due to the largeness of the configuration space, it is impractical to sample from it by direct computation of the probabilities. Markov chains Monte Carlo (MCMC) methods have been investigated by various researchers as an alternative to exact probability computation [1,2]. The general method is to simulate a Markov chain with the required probability distribution as its equilibrium distribution. If the chain is aperiodic and irreducible, the convergence is guaranteed.

The main goal is to introduce a method where a model can be defined and investigated. Based on this investigation, characteristics for this model (auto-Poisson) can be held and generalized for similar image patterns. Simulations of the process are achieved by using MCMC method, like Gibbs sampler, in simulated and real data.

2. Materials and Methods

Markov Random Fields Modeling

Spatial data under investigation can be defined as regular or irregular areas (Figure 1). Under these areas, spatial patterns cam be illustrated based on the structure of the images. The main goal is to model the spatial patterns in a way to represent as best as they can the real data starting from a simple model to a more complicated one (Figure 2).

Considering a 2-D space, which has been partitioned into n-pixels, labeled by the integers 1, 2, …n (rectangular space); x_ij or x_i denote the color for pixel (i) or (i,j); p(..|..) defined as the conditional probability distribution (Figure 3) ([3,4,5,6]).

If D is a finite lattice D = {(i, j), 1 ≤ i ≤ N, 1 ≤ j ≤ M}, then each pixel of the finite lattice D can be colored from the set {0, 1,…,c − 1}. Site j (¹ i) is said to be a neighbor of site i iff the functional form of p (x_i|x₁,…, x_i−1, x_i+1,…, x_n) is defined upon the variables x_j as $p\left(x_i|x_1,\dots ,x_{i-1},x_{i+1},\dots ,x_n\right)=p\left(x_i|x_{\partial i}\right)$, where ∂i is the set of pixels that are neighbous of pixel I, and x_∂i is the set of values of pixels that are neighbors of pixel i.

A neighborhood structure N = {N_i, ∀ i∈ S} is defined as a collection of subsets of S. The symmetry property is based on the following conditions: (i) i ∉ N_i (a site is not part of its neighborhood); (ii) j ∈ N_i ⇔ i ∈ N_j (i is in the neighborhood of j if and only if j is in the neighborhood of i). Define a nearest-neighborhood set as the set of sites with the property that p (x_ij|all other values) depends only upon the neighbors x_i−1,j, x_i+1,j, x_i,j−1, x_i,j+1 for each internal site (i,j) (Figure 4).

The first order model (four neighbors) is denoted by N = {(x_i,j, x_i−1,j), (x_i,j, x_i+1,j), (x_i,j, x_i,j−1), (x_i,j, x_i,j+1)} and second order (eight neighbors) is denoted by N = {(x_i,j, x_i−1,j), (x_i,j, x_i+1,j), (x_i,j, x_i,j−1), (x_i,j, x_i,j+1), (x_i,j, x_i−1,j−1), (x_i,j, x_i−1,j+1), (x_i,j, x_i+1,j−1), (x_i,j, x_i+1,j+1)} (Figure 5).

A clique is defined as a set that consists either of a single pixel or a collection of pixels that are neighbors of each other. Given a 2-D space of a subset S and neighbor structure ∂i, a clique is any set of pixels c ⊂ S, where for all i, j ∈ c, j ∈ ∂i. C can be defined as the set of all cliques. The concept of clique is directly combined with the calculation of the energy of the image, which can be considered as a statistical measure of the weight of the correlation between the pixels.

For the first order, the neighborhood structure is given by: {(i,j)}, {(i − 1,j), (i,j)}, {(i,j − 1), (i,j)}; and for the second order the neighborhood structure is given by: {(i,j)}, {(i − 1,j), (i,j)}, {(i,j − 1), (i,j)},{(i − 1,j − 1), (i,j)}, {(i − 1,j + 1), (i,j)}, {(i,j + 1), (i − 1,j + 1), (i,j)}, {(i,j + 1), (i − 1,j), (i,j)}, {(i − 1,j − 1), (i − 1,j), (i,j)}, {(i − 1,j + 1), (i − 1,j), (i,j)}, {(i − 1,j − 1), (i − 1,j + 1), (i + 1,j + 1), (i + 1,j − 1),(i,j)}.

A Markov random field (MRF) is a joint probability density on the set of all possible coloring values X on a finite lattice D, following the conditions ([3,4,7]): (i) p (x) > 0 for all x ∈ S. (Positivity). (ii) p (x_ij|all points (i,j)) = p (x_ij|neighbors). (Markovianity). (iii) p (x_ij|neighbors (i,j)) depends only on the configuration of the neighbors (Homogeneity). The probabilities on the condition (2) are called local characteristics. Condition (2) can be expressed as $p\left(x_i|x_j,i\ne j\right)=p\left(x_i|x_{\partial i}\right)$, where it is clear that we need a representative distribution (Gibbs). If p (x) is a jpdf of Χ_i ∈ S under any neighborhood model ∂i, the Gibbs distribution can be expressed based on the following form:

$$p\left(x\right)=\frac{1}{Z}\exp \left\{-\frac{1}{T}\left[U\left(x\right)\right]\right\}=\frac{1}{Z}\exp \left\{-\frac{1}{T}\left[{\displaystyle \sum _{c\in C}{V}_c\left(x_c\right)}\right]\right\}$$

(1)

where C are all the potential cliques, Ζ is the normalized constant/partitions function, Τ is is the temperature with Τ = 1, and U (x) is the energy function with $U\left(x\right)={\displaystyle \sum _{c\in C}{V}_c\left(x_c\right)}$ [8]. Under the Hammersley–Clifford theorem ([9,10]), if X is a discrete or continuous variable assigned to a pixel x representing a random field with neighborhood structure ∂i and jpdf p (x), then Χ is a MRF. iff p (x) can be expressed as a Gibbs distribution. The general form for the energy function can be defined by

$$U\left(x\right)={\displaystyle \sum _{i\in S}{V}_1\left(x_i\right)+{\displaystyle \sum _{i\in S}{\displaystyle \sum _{j\in S}{V}_2\left(x_i,x_j\right)}}}={\displaystyle \sum _{1\le i\le S}{x}_i{G}_i\left(x_i\right)+{\displaystyle \sum _{1\le i<}{\displaystyle \sum _{j\le S}{x}_i{x}_j{G}_{ij}\left(x_i,x_j\right)}}}+....$$

(2)

where G (.) is any arbitrary function ([3,4]). For example, for the second order neighborhood structure, the energy function is given by

$$U\left(x\right)={\displaystyle \sum _{i\in S}{V}_1\left(x_i\right)+{\displaystyle \sum _{i\in S}{\displaystyle \sum _{j\in S}{V}_2\left(x_i,x_j\right)}}}={\displaystyle \sum _{1\le i\le S}{x}_i{G}_i\left(x_i\right)+{\displaystyle \sum _{1\le i<}{\displaystyle \sum _{j\le S}{x}_i{x}_j{G}_{ij}\left(x_i,x_j\right)}}}$$

(3)

defined as pairwise interaction MRF models.

Specific spatial patterns can be described by particular models based on their neighbors, defined as auto-models. Assumptions for these models are: (1) The probability structure only depends on contributions of sites taken either as singular or in pairs. (2) The conditional probability distribution is a member of the regular exponential family of distributions $p\left(x_i|x_{\partial i}\right)=\exp \left\{A_i\left({\theta }_i\right)B_i\left(x_i\right)+C_i\left(x_i\right)+D_i\left({\theta }_i\right)\right\}$ where θ_i is a model parameter associated with site i and is a function of the values at sites neighboring site i. A_i (θ_i) can be defined as a potential interaction between pixels. Assuming the conditional probabilities and the pairwise only dependence between sites, the A_i (θ_i) must satisfy ([2,5,6]):

$$A_i\left({\theta }_i\right)={\alpha }_i+{\displaystyle \sum {\beta }_{ij}}{B}_j\left(x_j\right)$$

(4)

where β_ij = β_ji if i is neighbors with j and β_ij = 0 otherwise. As a final restriction, it is assumed that the function B_j (x_j) is linear in x_i with form $A_i\left({\theta }_i\right)={\alpha }_i+{\displaystyle \sum {\beta }_{ij}}{x}_j$. If β_ij = β_ji = β, it is an isotropic model, otherwise it is anisotropic. For the first-order system the models are: isotropic: A_i (α,β) = α + β (x_i−1,j + x_i+1,j + x_i,j−1 + x_i,j+1) and anisotropic: A_i (α, β₁,β₂) = α + β₁ (x_i−1,j + x_i+1,j) + β₂ (x_i,j−1 + x_i,j+1) (Figure 6a). For the second order, the models are: isotropic: A_i (α, β, γ) = α + β (x_i−1,j + x_i+1,j + x_i,j−1 + x_i,j+1)+γ (x_i−1,j−1 + x_i−1,j+1 + x_i+1,j−1 + x_i−1,j+1) and for anisotropic: A_i (α, β₁, β₂, γ₁, γ₂) = α + β₁ (x_i−1,j + x_i+1,j) + β₂ (x_i,j−1 + x_i,j+1) + γ₁ (x_i−1,j−1 + x_i−1,j+1) + γ₂ (x_i+1,j−1 + x_i−1,j+1) (Figure 6b) ([8]).

Considering the general notation of the auto-models, the potential interaction between pixels for the auto-Poisson model is given by

$$A_i\left({\lambda }_i\right)=\log \left({\lambda }_i\right)\Rightarrow {\lambda }_i=\exp \left(a_i+{\displaystyle \sum {\beta }_{ij}{x}_j}\right)$$

(5)

under the assumption that the conditional distribution of the pixels x_i given their neighbors has a Poisson distribution mean λ_i with form

$$p\left(x_i|x_{\partial i}\right)=\frac{{\lambda }_i^{x_i}\exp \left(-{\lambda }_i\right)}{x!}=\exp \left[\log \left({\lambda }_i\right)x_i-{\lambda }_i-\log \left(x_i!\right)\right]$$

(6)

The limitation of the model is that $\left|b\right|>\frac{1}{{\lambda }_i}$([2,5,6]), where b = mβc, c denotes number of colors, and m denotes the number of neighbors.

3. Results

Simulation Process Using MCMC Method

The investigation of the spatial patterns can be illustrated by the simulation procedures. A particular simulation process, where realizations from pseudo-samples can finally define the desired distribution, is the Monte Carlo Markov Chain (MCMC) ([11]). The goal of the process is to simulate the distribution p (x) using a particular realization X¹, X²,…, X^N on the Markov chain with transition probability. Under the process, asymptotic results can be archived where:

$$\begin{array}{l}{X}^t\underset{t\to \infty }{\overset{d}{\to }}X~p(x) ; \frac{1}{t}{\displaystyle \sum _{i=1}^{t}f(x^i)\underset{t\to \infty }{\to }}{E}_p\left\{f(x)\right\}\\ \end{array}$$

(7)

with the expectation E_p {f (x)} under estimation. The corresponding empirical average will be given by

$${\overline{f}}_N=\frac{1}{N}{\displaystyle \sum _{t+1}^{N}f\left(x^{(t)}\right)}$$

(8)

Gibbs sampler is a special case of MCMC methods. Consider Χ_i ∈ S a discrete or continuous value for a random field in a rectangular lattice system S with neighborhood structure ∂i. If X is defined as a pdf p (x) based on Gibbs distribution, then the conditional probability of a value given its neighbors can be defined as

$$p(x_s|x_j,s\ne j)=\frac{\exp \left\{-U\left(x_s|x_j,s\ne j\right)\right\}}{{\displaystyle \sum _{x_s{}^{*}\in S}\exp \left\{-U\left(x_s^{*}|x_j,s\ne j\right)\right\}}}$$

(9)

Each value could be replaced by the conditional probability

$$p(x_s|x_j^k,j\ne s)=p(x_s|x_{\partial s})$$

(10)

The algorithm stopped when all the pixels were replaced. The Gibbs sampler algorithm is given by the following pseudo-algorithm ([12,13,14,15]):

Step 1: Chose randomly pixel x_i

Step 2: For t = 0 until ∞

Replace the value x^(t) with x^(t+1) based on

$$X^t\underset{t\to \infty }{\overset{d}{\to }}X~p(x) ; \frac{1}{t}{\displaystyle \sum _{i=1}^{t}f(x^i)\underset{t\to \infty }{\to }}{E}_p\left\{f(x)\right\}$$

(11)

Step 3. Continue the process until all the pixels have been replaced

Realizations for the auto-Poisson model are given in Figure 7 with parameter a. First order isotropic with α = 0.75 and β = −0.25; b. Second order isotropic with a = −0.93, β = −0.33, and γ = −0.37. Finally a comparison between biological images from microscopes with realizations from the simulated auto-Poisson model considering the first-order isotropic model (α = −1, β = −2) is given in Figure 8, where it is clear that both images have similar patterns.

4. Conclusions

Pattern analysis is a process where regular synthesis or patterns can be recognized. Most of the time, these patterns can be analyzed using texture models based on the spatial structure of the images. The spatiality of the images can be defined by considering the homogeneity of the regions. These regions might be represented by considering identity patterns under the neighborhood structure of the image. Models for the estimation of these spatial patterns could be introduced considering the local characteristics of the image, leading us to Markov random fields (MRF) models. These patterns could be used to simulate real phenomena like biology, ecology, or medicine. In this work, a presentation of Markov random field models (MRF) considering specific conditional distributions like auto-models was analyzed. A specific auto-model under Poisson distribution (auto-Poisson) were defined and simulations using MCMC methods as Gibbs samples were illustrated using simulated and real images patterns. Modification of auto-Poisson models was successfully performed in many case like social networks [16], spreading disease modeling [17,18], ecological models [19,20,21,22]. In these works, extensions of Markov random fields (MRF) models were applied using network analysis considering the neighborhood structure of the pixels. Under this extension, homogeneous regions could be illustrated considering the connections between similar pixels. Another important application is based on the mapping analysis where simulations of auto-Poisson models might be used to reconstruct estimated maps under various regions, analyzing phenomena like cancer diseases [23]. Last but not least, the proposed models could be used for firm location analysis [24] to estimate the optimal positions for firm’s establishment.