A study on application of Markov random fields in image restoration and its efficiency

Image restoration has been a popular and active field of research for decades. Images are destroyed when exposed to ‘noise’ which can occur due to physical contact or electrical/electronic interference. Here, Bayesian statistical techniques and Markov random field (MRF) theory were used to restore a black and white binary image corrupted by additive Gaussian noise with zero mean and constant variance. The binary image was used as a Markov random field. An image is comprised of pixels and these pixels have a regular two-dimensional lattice structure. A matrix including ±1 values was generated randomly. It was represented and interpreted as an Ising model in Statistical Mechanics. The probability distribution of the Ising model was used as the prior distribution (a Gibbs or Boltzmann distribution). Likelihood function was obtained by using the random matrix and the observed corrupted image. Markov Chain Monte Carlo (MCMC) method was used to simulate posterior distribution which again turns out to be a Gibbs or Boltzmann distribution. More specifically, Metropolis-Hastings algorithm which is one of the popular MCMC algorithms was used in this simulation. In this study, Peak Signal to Noise Ratio (PSNR) and Mean Squared Error (MSE) methods were used to measure the quality of restored images. MATLAB (R2013a (8.1.0.604)) was used to construct the program in this research work. Finally, high quality images were restored which were almost similar to the original image. A decrease in restored image quality was observed with the increase in noise. When the image size was increased, a higher number of iterations was required to obtain an acceptable level of quality in the restored image.


Introduction
Image restoration has been proven to be an active field of research today. The development of effective methods of denoising an image still remains a big challenge (M. Kaur & Behal, 2013). Images are destroyed due to various physical phenomena. On a mathematical context, destroyed images have 'noise' which can occur due to physical contact or electrical/electronic interference. The purpose of image denoising is to maintain the main features of the original image as much as possible and remove noise from the corrupted image (Saxena et al., 2017). Image processing has been very useful in computer vision problems, medical imaging techniques, satellites, radar and sonar, etc. There are several noise types that affect image corruption such as image noise as impulse noise (salt-and-pepper noise), Gaussian noise, uniform noise, Brownian noise (fractal noise), Rayleigh noise, Poisson noise, etc. (A. Kaur & Chopra, 2012;Raval & Gagnani, 2012). In this research, Gaussian noise with zero mean and constant variance was used to obtain a corrupted binary (Black and White) image.
Restoring original image from the corrupted image has been studied by a number of authors using different techniques. Over the last twenty years, several techniques have been developed to support various binary image processing applications. According to Open access article RESEARCH ARTICLE Journal of the University of Ruhuna, Sri Lanka 8(1), 2020 previous literature, novel minimization model (J. Zhang & Ye, 2011), positive semidefinite programming (Shen et al., 2006), signomial programming optimization approach (Shen et al., 2007), neural network models (Su et al., 1993) and discrete tomography (Nemeth & Balazs, 2013) models have been used in binary image restoration. Here, Bayesian statistical techniques and Markov random field theory were used. Markov random field method based Ising model and 1st order neighborhood system with Markov property were used to generate a prior distribution. Ising model can be mapped to the binary image. The probability of a particular image pixel depends on the probabilities of the neighboring pixels. Therefore, the consideration of neighborhood system is significantly important in image restoration. Bayesian technique gives posterior distribution of the restored image which is very similar to the real image. In this study, images were restored in two ways that is by increasing noise variance ( ) while keeping the image size constant and by increasing the image size while keeping noise variance ( ) constant.

Grayscale Image to Binary Image
A binary image was used as the image source. White image pixel values were set to +1 and black image pixel values were set to −1 in computation by using MATLAB software package. A grayscale image has its pixel intensity of black information from 0 to 1. Therefore, image pixels were compared with 0.5 value. If pixel value was above 0.5, then it was set to +1 and if pixel value was below 0.5, then it was set to −1. That image was used as the binary image.

Introduction to Markov Random Fields (MRFs)
The Markov random field is an -dimensional random process defined on a regular lattice structure. Generally, the lattice structure is a regular 2dimensional (2-D) grid in the plane. Suppose = { 1×1, . . . , × } is a family of random variables defined on the state space and f = 1×1 , . . . , × is called a configuration of Markov random field . Here, × is the dimension of the 2-D grid plane. × = is also called the set of all sites.
are the neighbors of the ℎ site.

Gibbs Random Field
A collection of random variables is called a Gibbs random field on state space with respect to neighborhood system if and only if their configuration follow a Gibbs distribution. A Gibbs distribution can be written as is the normalizing constant. It is also called the partition function. Here, , and ( ) are temperature, Boltzmann constant and energy function respectively (Li, 2009).

Ising Model
Ising model is the simplest classical model for describing the magnetic behavior of atoms. It can also be interpreted as a MRF. It is a system of interacting particles, arranged in a regular lattice structure. Particles have two magnetic spin orientations that is spin up (+1) and spin down (-1). Adjacent particles that have the same spin (-1, -1) or (+1, +1) are in a lower energy state than those with opposite spins (+1, -1) or (-1, +1). Total energy of the system depends on the orientation of a particular spin compared to its neighbors. Given the spin orientations of all particles in the system, one may compute the total energy. If the variable denotes the spin of particle , then the total energy of the system (Hamiltonian of the system) for each configuration is where and are parameters which correspond to the energy associated with nearest-neighbor interactions and interactions with the external field, respectively. The parameter corresponds to the presence of an "external magnetic field".
Suppose that the set represents all possible configurations, and ( ) is the energy of configuration . Statistically, the probability density function of any particular configuration is described by using a Gibbs distribution (Cipra, 1987) of the form We can use the model (6) to describe the interaction between black and white pixels in binary image. Figure 2 illustrates disagreeing edges. The edge between two pixels with different values is a disagreeing edge which is colored as bold red line here. In this figure, number of disagreeing edges are 8. Disagreeing edges have higher energy than other edges. Therefore, we obtain the density function of configuration in image as, where = 1 , is the number of disagreeing edges and is a constant parameter. can be estimated by binary black and white image pixels. In literature, the value of has been taken as 1.0266 (Xiong & Hong, 2017).

Introduction to Bayesian Statistics
In real situations, we have noisy images only. We have to restore original image from the noisy image. The probabilities of the neighboring pixels affect the probability of a particular image pixel. It can be used as the probability model. Bayesian statistical technique can be used to estimate the original image by its posterior distribution. Therefore, Bayesian statistical technique is an efficient simulation method which allows the restoration of noisy images.
A sampling model for the observed data conditioned on the unknown parameter such that ∼ ( | ) is introduced where ( | ) is the likelihood function for observed data.
Initially, the prior distribution ( ), such that ∼ ( ), is specified. The prior distribution ( ) gives our belief about the possible values of the parameter before taking the data. Bayes' theorem can be
The objective is to estimate the true classification while the classification is given. Metropolis-Hasting algorithm is used to simulate the posterior of true classification given . It is denoted by ( |Y). Metropolis-Hasting algorithm is one of the Markov Chain Monte Carlo (MCMC) algorithms. Finally, a posterior distribution ( | ) that converges to the true posterior distribution of ( | ), can be simulated by using Metropolis-Hasting algorithm.

Metropolis-Hasting Algorithm for Ising Image Model
We begin with the space of all configurations in which each configuration is denoted by the matrix with the indexing ( , ) ⟼ = ( − 1) × + , the Metropolis-Hasting algorithm would include the following steps:  (Xiong & Hong, 2017).
The likelihood function of model is given by Then likelihood ratio is obtained as Journal of the University of Ruhuna, Sri Lanka 8(1), 2020 where ′ = − .
Priors' ratio can be obtained using Ising model distribution.
where = Constant also called interaction strength of the image.
= Number of disagreeing edges on ′ = Number of disagreeing edges on ′

Image quality metrics
Image quality metrics are used to measure the quality of image after a restoration process comparing the original and restored image. Measuring image quality is essential for many image processing applications. Peak signal to noise ratio (PSNR) is the most popular and widely used technique of the image quality metric. Further, mean square error (MSE) technique is also used.

Mean Square Error (MSE)
The MSE is the cumulative squared error between the compressed and the original image. Small value of MSE improves image quality and reduces the error.

Peak Signal to Noise Ratio (PSNR)
PSNR is used to measure the quality of any reconstructed, restored or corrupted image with respect to its reference or original image. PSNR is measured in decibels ( ). The higher value of PSNR indicates higher quality of image (Silpa & Mastani, 2012). PSNR is obtained by comparing the mean square error (MSE) with the maximum possible value of the pixel (L) as follows: where is the original image, is the denoised or restored image and is the dimension of the image. (for gray scale image (8 bit per pixel) = 2 8 − 1 = 255 and binary image (1 bit per pixel) = 2 1 − 1 = 1.)

Results and Discussion
Restored images were obtained in two different ways: Increasing noise variance ( ) while keeping a constant image size and increasing the image size while keeping a constant noise variance.

Increasing noise variance ( ) while keeping a constant image size ( × )
The following Figure 3 shows the obtained restored images while is increased from 0.5 to 2.5 by 0.5 intervals. It was obvious that with the increase in noise the restored image quality decreased. When noisy images with variance = 0.5, and = 1.   image are 0.9950 and 0.1000 respectively and the MSE reduction percentage value of the restored images is 89.95%.
The Figure 6 shows the similarity between noisy and restored images with original image (image size = 1048 × 750). MSE values of noisy image and restored image are 1.0003 and 0.0426 respectively while the MSE reduction percentage value of the restored images is 95.74%.
The Figure 7 shows the similarity between noisy and restored images with original image (image size = 2096 × 1500). MSE value of noisy image is 1.0005 and that of restored image is 0.0227. MSE reduction percentage value of the restored images is 97.73%. Figure 8 shows the behavior of the image sizes and PSNR values. At 5000 × 10 4 iterations, restoring images converged. When the image size is decreased, the restored converged image can be obtained rapidly but with a lower quality. The observed PSNR values of three different sizes of noisy ( = 1) images 524 × 375, 1048 × 750, and 2096 × 1500 pixels are 9.9996 , 13.709 and 16.4434 respectively.

Conclusion
According to the results, when noise ( ) is increased from 0.5 to 2.5 by 0.5 intervals, the converged restored image quality decreases as MSE reduction percentage values of the restored images are 97.80%, 95.58%, 93.14%, 89.27% and 88.15% respectively. It can be concluded that, the quality of the restored image is very high and the observed images are very similar to the original image. When noise variance ( ) is increased keeping a constant image size, the MSE values relatively increase while the PSNR values relatively decrease.
When increasing the image size (524 × 375, 1048 × 750, and 2096 × 1500) keeping a constant noise variance ( = 1), the converged restored image quality increases as MSE reduction percentage values of the restored images are 89.95%, 95.74% and 97.73% respectively. It can be concluded that, when the image size is increased, converged restored image quality increases. Furthermore, the number of iterations required for obtaining the converged restored image also increases.
In this study, 1 st order neighborhood system was used to get Markov property in the Ising model. For the particular study, 2 nd or higher order neighborhood systems can be used given that it provides more accurate and better quality results for the restored images.