Unsupervised segmentation of images using bi-dimensional pairwise Markov chains model

1.   Introduction

The hidden Markov chain (HMC) model, as characterized in ^[1,2,3,4,5], is a doubly stochastic process, usually denoted by $(\textbf{X}, \textbf{Y})$. The process $\textbf{X}$, which models a to be segmented image, is assumed to be hidden and Markovian, and the process $\textbf{Y}$, which models the noisy version of the image to be recuperated, is assumed to be observed and real. The success of the HMC model stems from the fact that, when the noise is not too complex, several Bayesian classification approaches such as maximum A Posteriori (MAP) or maximal A Posteriori (MPM) estimators can be used to determine the image $\textbf{X} = \textbf{x}$ from the observed version $\textbf{Y} = \textbf{y}$. In addition, the spatial regularity of pixels may hinder the effectiveness of the HMC model in image segmentation.

We consider the problem of segmenting a satellite image into two classes, "water" and "trees". The hidden image represents all of the classes, whereas the observed image represents all of the gray levels on each pixel. In general, all neighboring pixels are considered to have a higher probability of the same class than pixels situated further from each other. Nevertheless, in this situation, we can detect pixels that are among a set of pixels in identical groups or near edges but have distinct appearances. According to the classical hypothesis, this attribute cannot be accommodated in the HMC model. For this reason, the HMC model has been augmented to include the pair-wise Markov chain (PMC) model The PMC model is more general than the HMC model, given the process $\textbf{X}$ is not necessarily a Markov chain and the pair $(\textbf{X}, \textbf{Y})$ is directly Markovian. Another advantage of PMC against HMC is its capacity to easily take into account the noise correlation. Furthermore, the PMC model, like HMC, is capable of implementing Bayesian MAP and MPM restorations.

Another major problem is that when a one-dimensional sequence is generated from a two-dimensional image using any reading process, such as human preference score (HPS), some locality information of the pixels is lost. What prompted us in this study was the success of adding additional components to the observations in the studies ^[6,7,8] in the framework of the HMC model by taking advantage of a neighborhood of each pixel (for example, four nearest neighbors), which culminated in an impressive model on a range of synthetic images. In this work, we apply the same concept to the PMC model in order to tackle the main issue of spatial regularity of pixels, as well as the problem of pixel information loss in an image, which is the objective of the segmentation process. As a result, we formulate a PMC model with bidimensional observed data (BPMC), whose most important task is to segment images while taking into account the above-mentioned last two difficulties. The focus of this study is to look into the potential benefits of BPMC over PMC in terms of unsupervised segmentation. Initially, the unidentified parameters of the model must be estimated. In the case of PMC, we cannot simply maximize the likelihood analytically, which is the case for the expectation-maximization (EM) algorithm ^[9,10], at each iteration, so we need to apply iterative conditional estimation (ICE) and SEM ^[11,12,13]. In general, the ICE and structural equation models (SEM) algorithms are very powerful when used with Gaussian estimation and generalized mixture. In this study, we adapt the ICE and SEM algorithms to the new model. Furthermore, SEM and ICE algorithms are examined in the presence of Gaussian noise for both PMC and BPMC, and it is discovered that SEM and ICE algorithms are more adaptable because likelihood is not always required. In addition, the estimators beyond the maximum likelihood can also be utilized. Therefore, we need to answer the following questions.

● Does the use of BPMC significantly improve restoration results, compared to the PMC model, and how do ICE and SEM algorithms work in the BPMC context?

● Does the use of BPMC compared to PMC permit the image to be converted into a one-dimensional sequence without losing pixel locality information?

The last point ought to be considered since all kinds of hidden Markov chain models based on HPS have shown potential for image segmentation, and in certain situations, they may even be competing with the hidden Markov field (HMF) model in classification ^[14,15,16]. On top of that, the HMF model allows for a more precise and more intuitive modeling of the spatial relationships between pixels ^[17], while the PMC model requires HPS to account for the spatial details of the pixels. Conversely, the BPMC model takes into account temporal and spatial pixel information by incorporating an additional element right into the observed process, allowing it to compete with the PMC model in this regard.

The rest of this paper is structured into seven sections. Section 2 presents studies relevant to the current study. In Section 3, we provide the reader with a mathematical background of the BPMC model. Section 4 presents the two models discussed in this study and their properties. Section 5 is committed to the simulation of PMC and BPMC models in a supervised way using MPM classifier-based classifications. Section 6 deals with estimating parameter algorithms; here, we describe the main steps of the ICE and SEM algorithms in the case of the two models, PMC and BPMC. In Section 7, with numerical experiments for various noise factors, we give a study of comparison of the two estimation procedures based on the accuracy of parameter estimates for both models. Furthermore, using error rates and the peak signal-to-noise ratio (PSNR) index, we demonstrate the performance of these models on a sample of synthetic images. Finally, conclusions and some perspectives are discussed in Section 8.

2.   Literature review

The pair-wise terminology was first attached to the hidden Markov model in the works ^[18,19,20], where the main motivation was to solve the problem of correlation between noise and pixels located on the image boundary. In fact, this difficulty has been emphasized for the first time in the pair-wise Markov field (PMF) model and the pair-wise Markov tree (PMT) model. The PMF model has been successfully used for textured image segmentation ^[20], real synthetic-aperture radar (SAR) images ^[21], and hyper-spectral images ^[22]. The PMT model has also been used for image and signal restoration tasks ^[19]. Recently, the PMC model ^[18] has been used for applications including signal and image processing. The PMC model has specific applications in image processing, such as textured and real images segmentation ^{[14,23,24,25,26,27]} and text segmentation ^[28,29], where the authors used an HPS for converting the image and copula to model conditional densities of class. Images with multi-scale can be segmented using PMC ^[30]. Other uses include fuzzy PMC to treat spatially correlated fuzzy classes between observations ^[31,32].

Several other methods based on hidden Markov models have been attempted in order to solve the question of the relationship between pixels in the chain after the transformation of a bidimensional image. Among these, an approach suggested by ^[15] is based on a combination of the HMC and HMF models. The first model is employed in the estimation phase and the second is exploited in the final classification to determine the $\textbf{X}$ configuration. Besides, another strategy in ^[33,34] is centered on the integration of spectral and contextual information into the HMC model, where the structure of the gathered data is no more Markov when modeling an image using HPS. Additionally, the study ^[35] entails another suggestion for an algorithm that incorporates fuzzy C-means combined with HPS. In ^[36], a new HMC is proposed for representing the semantic and spatial interactions between pixels.

A new segmentation method in ^[7] has recently been introduced that tackles the issues of spatial and temporal pixels via contextual HPS. This method adds an additional element to the observation process in the HMC model so that each pixel recovered from the HPS is linked to another pair of pixels next to it in the image but not in the chain. Similarly, in ^[8], a second component is introduced in the process Y and a two-dimensional model is generated with estimation and segmentation techniques adapted to this scenario. Table 1 summarizes all of the methods covered here.

3.   Mathematical context

In this part, we describe a brief mathematical basis for the BPMC model given in this study. We consider a couple of stochastic processes $(\textbf{X}, \mathbf{Y}) = ((X_n)_{1\leq n\leq N}, (Y_n)_{1\leq n\leq N})$, where $N$ is the number of pixels. As previously explained, in the segmentation of images with the hidden Markov model (HMM), we need to estimate the latent process $\textbf{X} = \textbf{x}$ via the observed process $\textbf{Y} = \textbf{y}$.

In the context of the PMC model with two-dimensional observed data, we consider another realization of a random variable $\bar{\textbf{Y}} = \bar{\textbf{y}}$ modeling information about each pixel based on its neighborhood in the image, but not in the chain, in such a way that the observed process is two-dimensional and is denoted by $(\mathbf{Y}, \bar{\textbf{Y}})$. Then, we intend to estimate the parameters of the joint probability $p(\textbf{x}, \textbf{y}, \bar{\textbf{y}})$ by using data unsupervisedly. Based on the observation, Bayesian estimators are used to restructure the hidden process. In the present model's segmentation problem, it considers two laws: An a priori law $p(\textbf{x})$ and a two-dimensional conditional density $p(\mathbf{y}, \bar{\textbf{y}}|\textbf{x})$.

To take on the above-mentioned segmentation problem, the a posteriori distribution must be identified, which contains all of the information on $\textbf{X}$ that is available at $(\mathbf{Y}, \bar{\textbf{Y}})$. This distribution was supplied by: $p(\textbf{x}|\textbf{y}, \bar{\textbf{y}}) = \frac{p(\textbf{x}, \textbf{y}, \bar{\textbf{y}})}{ \sum\limits_{\textbf{x}}p(\textbf{x}, \textbf{y}, \bar{\textbf{y}})}$, where $p(\textbf{x}, \textbf{y}, \bar{\textbf{y}}) = p(\textbf{x})\times p(\textbf{y}, \bar{\textbf{y}}|\textbf{x})$.

Let $L$ be a bivariate loss function and $\hat{\textbf{x}}$ be an approximation of the hidden realization, where the latter is obtained by this method:

In Bayesian statistics, the two best-known estimators are MAP and MPM, which are associated respectively with the two loss functions $L_{MAP} (\textbf{x}, \textbf{x}') = L(\textbf{x}, \textbf{x}')$ and $L_{MPM}(\textbf{x}, \textbf{x}') = \sum\limits_{n = 1}^{N}L(x_n, x'_n)$, with $L$ here representing the Kronecker symbol.

The two estimators are consequently and respectively defined by:

and

When comparing these two estimators, MPM maximizes the a posteriori marginal probability for each $x_n$, while MAP maximizes the a posteriori probability $p(\textbf{x}|\textbf{y}, \bar{\textbf{y}})$ in a direct way, and allows us to estimate the whole $\textbf{x}$ sequence. Next, we only utilize the MPM estimator since it is straightforward to use. We use iterative strategies to estimate the joint distribution of the unobserved process using the data available $(\textbf{y}, \bar{\textbf{y}})$. For calculating the joint distribution, we often use a conditional expectation approximation in the case of the ICE method or stochastic generators to estimate a sequence of model parameters from the dataset in the case of the SEM algorithm.

4.   Models

Let's consider two stochastic processes $\textbf{X} = (X_n)_{1\leq n\leq N}$ and $\textbf{Y} = (Y_n)_{1\leq n\leq N}$, where $N$ is the number of pixels of an image. The process $\textbf{X}$ is the unknown image whereas $\textbf{Y}$ is the observed one. Each random variable $X_n$ has the values from a finite set of $K$ classes $\Omega = \{\omega_1, \ldots, \omega_K\}$, and each $Y_n$ has attributed the values from the set of real numbers $I\!\!R$. Realizations of $\textbf{X}$ and $\textbf{Y}$ are denoted by $\textbf{x} = (x_n)_{1\leq n\leq N}$ and $\textbf{y} = (y_n)_{1\leq n\leq N}$, respectively.

We suppose that the process $\textbf{X}$ is stationary and that the random variables $\textbf{Y} = (Y_n)_{1\leq n\leq N}$ are correlated and conditionally depend on $\textbf{X}$. As in all Markov models, the difficulty of classification is estimating the latent process $\textbf{X} = \textbf{x}$ based on the observed process $\textbf{Y} = \textbf{y}$.

4.1.   PMC model

This study uses the PMC model and proposes updating the classical pair-wise model. The classical PMC model is a doubly stochastic process, usually noted ($X$, $Y$), where $Y$ (one-dimensional) represents the pixels of a hidden image $X$. The same concept for the bidimensional PMC model is true, except for the process $Y$. In the classical model, a direct Hilbert Peano scan is used to transform the image from its two-dimensional version to a one-dimensional chain, followed by an inverse scan to reconstruct the final image, but in the suggested model, the second component of $Y$ provides for this goal. The visual results revealed a significant advantage of the proposed model over the classical model, especially for isolated pixels in the to-be-segmented image.

An image of a bidimensional form can be transformed into a one-dimensional form inevitably by combining the image column by column or line by line. Moreover, we can use the classical HPS for an image $2^p\times2^p$ ^[2], or the generalized HPS for any given image ^[37]. In a moderately different sense than above, now let's consider the pair-wise process $\textbf{Z} = (Z_1, \ldots, Z_n, \ldots, Z_N)$ corresponding to $\textbf{X}$ the hidden process and $\textbf{Y}$ the observed process, where $Z_n = (X_n, Y_n)$ for $n = 1, \ldots, N$. Lowercase letters are used to identify the realizations of such variables and processes. To keep things straightforward and brief, we use the letters $p$ and $f$, respectively, to describe various distributions on $\Omega$ and densities on $I\!\!R$ in the following text. $p(z_n)$ is also employed for $p(Z_n = z_n)$. Referring to the previous, the pair $\textbf{Z}$ is a PMC model when its joint distribution can be given as follows $p(\textbf{z}) = p(z_1) \prod\limits_{n = 2}^{N}p(z_n|z_{n-1})$, with

and

Based upon whether the process $\textbf{X}$ is stationary, we can also define the $\textbf{Z}$ distribution only as $p(z_1) p(z_2|z_1) = p(z_1, z_2)$. Thus the following proposition is inspired by the following works ^[4,6,38].

Proposition 4.1. For $\boldsymbol{Z}$ a stationary PMC model related to processes $\boldsymbol{X}$ and $\boldsymbol{Y}$, and for an $n$ from $\{2, \ldots, N\}$, the following results are obtained.

(1) Both the distributions of $\boldsymbol{X}$ conditionally to $\boldsymbol{Y} = \boldsymbol{y}$ and of $\boldsymbol{Y}$ conditionally to $\boldsymbol{X} = \boldsymbol{x}$ are Markovians.

(2) If $p(x_{n}|x_{n-1}, y_{n-1}) = p(x_n|x_{n-1})$ and $p(y_n|x_n, x_{n-1}, y_{n-1}) = p(y_n|x_n, y_{n-1})$, then $\boldsymbol{Z}$ is an HMC with dependent noise.

(3) If $p(x_{n}|x_{n-1}, y_{n-1}) = p(x_n|x_{n-1})$ and $p(y_n|x_n, x_{n-1}, y_{n-1}) = p(y_n|x_n)$, then $\boldsymbol{Z}$ is an HMC with independent noise.

(4) If $p(y_{n-1}|x_{n-1}, x_n) = p(y_{n-1}|x_{n-1})$, then $\boldsymbol{Z}$ is an HMC with dependent noise.

Proof. (1) Let $\textbf{Z} = (\textbf{X}, \textbf{Y})$ be a PMC model. Then, $n$ belongs to $\{2, \ldots, N\}$, and we have

We also have

This gives $p(x_{n}|x_1, \ldots, x_{n-1}, \textbf{y}) = p(x_{n}|x_{n-1}, \textbf{y})$, therefore $\textbf{X}$ conditionally to $\textbf{Y} = \textbf{y}$ is a Markov chain. In the same way, we can show that $\textbf{Y}$ conditionally to $\textbf{X} = \textbf{x}$ is a Markov chain, and $\textbf{x}$ and $\textbf{y}$ do a symmetrical function here.

(2) $\textbf{Z}$ is a PMC model, then

from which $\textbf{Z}$ is an HMC with dependent noise, identified by the initial distribution $p(z_1)$ by

and the transition matrix $p(z_n|z_{n-1})$ given by

(3) $\textbf{Z}$ is a PMC model, then

From which $\textbf{Z}$ is an HMC with independent noise.

(4) Let us show here that the process $\textbf{X}$ is a Markov chain. We have, for $n$ from $\{2, \ldots, N\}$,

Taking advantage that $p(y_{n-1}|x_{n-1}, x_n) = p(y_{n-1}|x_{n-1})$, we thus have $p(x_n|x_1, \ldots, x_{n-1}) = p(x_n|x_{n-1})$, which completes the proof.□

Remark 4.1. Properties $2$ and $3$ in ${Proposition \; 4.1}$ demonstrate the strict generality of the PMC model compared with the HMC model.

4.2.   BPMC model

In this model, we consider another observed process noted $\bar{\textbf{Y}}$ given by $\bar{\textbf{Y}} = (Y_1, \bar{Y}_1, \ldots, Y_n, \bar{Y}_n, \ldots, Y_N, \bar{Y}_N)$. Consequently, for each $n = 1, \ldots, N$, the random variable $X_n$ is associated with the pair $(Y_n, \bar{Y}_n)$, where $\bar{Y}_n$ is related only to the observation $Y_n$ as defined below. The problem remains unchanged; determining $\textbf{X}$ based on $\textbf{Y}$. To establish this new model, we may then state the following conditions (Figure 1).

(1) The process $\textbf{Z} = (\textbf{X}, \textbf{Y})$ is a PMC model.

(2) The random variables $(\bar{Y}_n)_{1\leq n\leq N}$ are conditionally independent to $\textbf{Z}$.

(3) Each $\bar{Y}_n$ is conditionally distributed to $\textbf{Z}$ and its distribution is conditionally to $Z_n$.

The model to be studied $(\textbf{Z}, \bar{\textbf{Y}})$, called the BPMC model, could be interpreted as a semi-hidden Markov chain model ^[4]. Its distribution can be stated as follows.

Let's emphasize that the BPMC model is particularly interesting because, given the pair $(\textbf{Y}, \bar{\textbf{Y}}) = (\textbf{y}, \bar{\textbf{y}})$, the distribution of $X$ is a Markov. Using the same logic as the demonstration of the first point in ${Proposition \; 4.1}$, we get the following proposition.

Proposition 4.2. For $(\boldsymbol{Z}, \bar{\boldsymbol{Y}})$, a stationary BPMC model related to processes $\boldsymbol{X}$ and $\boldsymbol{Y}$, the following result is obtained.

The distributions of $\boldsymbol{X}$ conditional on $(\boldsymbol{Y}, \bar{\boldsymbol{Y}}) = (\boldsymbol{y}, \bar{\boldsymbol{y}})$ is a Markov chain.

5.   Restoration problem of simulated models

The reason behind the achievements of the HMC model can be attributed to the fact that the process $\textbf{X}$ conditionally on $\textbf{Y} = \textbf{y}$ is a Markov chain and its transitions can be computed. This observation holds true in the case of the PMC model. However, when it comes to the BPMC model, this conditional distribution exhibits a complex structure. That is precisely why we employ the distribution $p(\textbf{x}|\textbf{y}, \bar{\textbf{y}})$. Both models, PMC and BPMC, are assumed to be stationary, such that their distributions are given by $p(x_1, y_1, x_2, y_2)$ for PMC and $p(x_1, y_1, \bar{y}_1, x_2, y_2, \bar{y}_2)$ for BPMC.

In this section, we contrast the two models using the MPM restoration algorithm in the Gaussian context. Both models utilize Gaussian densities, which implies that $p(y_{n-1}, y_n|x_{n-1}, x_n)$ and $p(\bar{y}_{n-1}|x_{n-1})$, for each $n$, are Gaussian (in the BPMC model case, we choose $p(\bar{y}_{n-1}|x_{n-1})$ to be equal to $p(\bar{y}_{n-1}|x_{n-1}, x_n)$). In this scenario, we have $K^2$ probabilities $p(x_{n-1} = \omega_i, x_n = \omega_j)$, $K^2$ mean vectors $\mu_{ij} = \binom{\mu^1_{ij}}{\mu^2_{ij}}$, and $K^2$ variance-covariance matrix $\Sigma_{ij} = \left(\begin{array}{cc} \sigma^{1}_{ij} & \sigma^{12}_{ij} \\ \sigma^{12}_{ij} & \sigma^{2}_{ij} \\ \end{array}\right)$ of the bidimensional densities $f_{ij}(y_{n-1}, y_n)$ to represent the PMC model. Additionally, the BPMC model is represented by the parameters $\mu^1_{ij}$, $\mu^2_{ij}$, $\sigma^{1}_{ij}$, and $\sigma^{2}_{ij}$ of the $K^2$ mono-dimensional densities $f_{ij}(\bar{y}_{n-1})$. In the following subsection, we will simulate the two processes $\textbf{X}$ and $\textbf{Y}$ assuming that all the aforementioned parameters are initially known. Furthermore, we will reconstruct the process $\textbf{X}$ only based on the process $\textbf{Y}$ using the PMC model, while considering both $\textbf{Y}$ and $\bar{\textbf{Y}}$ simultaneously using the BPMC model.

5.1.   Simulation of processes in models

Here, we evoke the previously defined distributions of the PMC model and BPMC model. The distribution of the PMC model is given by $p(x_1, y_1) = p(x_1)p(y_1|x_1)$ and $p(x_n, y_n|x_{x-1}, y_{n-1}) = p(y_n|x_{n-1}, y_{n-1}, x_n) p(x_n|x_{n-1}, y_{n-1})$, while the distribution of the BPMC model is represented by $p(x_1, y_1, \bar{y}_1) = p(x_1)p(y_1|x_1)p(\bar{y}_1|x_1, y_1)$, $p(x_{n}, y_{n}|x_{n-1}, y_{n-1})$, and $p(\bar{y}_{n}|x_{n}, y_{n})$.

We can simulate the processes $\textbf{X}$, $\textbf{Y}$ and $\bar{\textbf{Y}}$ of models concurrently and interchangeably using the following approach.

● For the PMC model, $x_1$, $y_1$, $x_n$, and $y_n$, for each $n = 2, \ldots, N$, are simulated with the following probabilities:

● For the BPMC model, $x_1$, $y_1$, $x_n$, and $y_n$ are simulated in the same way as above. Additionally, we have simulated $\bar{y}_{1}$ and $\bar{y}_{n}$, for each $n = 2, \ldots, N$, with these probabilities:

Remark 5.1. It can be seen here that (to dig deeper, see ^[14,39])

● $f_{x_1, x_2}(y_1)$ is a mono-dimensional Gaussian density $\mathcal{N}\left(\mu^1_{ij}, \sqrt{\sigma^1_{ij}}\right)$.

● $f_{x_{n-1}, x_{n}}(y_{n-1})$ is a mono-dimensional Gaussian density $\mathcal{N}\left(\mu^1_{ij}, \sqrt{\sigma^1_{ij}}\right)$.

● it can be shown that $y_n$ is generated according to a Gaussian drawing, where $p(y_n|x_n, x_{n-1}, y_{n-1})$ is a mono-dimensional Gaussian density of mean $\mu^2_{ij}+\frac{\sigma^{12}_{ij}}{\sigma^1_{ij}}(y_{n-1}-\mu^1_{ij})$ and standard deviation $\sqrt{\frac{\sigma^1_{ij}\sigma^2_{ij}-\sigma^{12}_{ij}}{\sigma^1_{ij}}}$.

● $f_{x_1, x_2}(\bar{y}_1)$ is a mono-dimensional Gaussian density $\mathcal{N}\left(\mu^1_{ij}, \sqrt{\sigma^1_{ij}}\right)$.

● $f_{x_{n-1}, x_{n}}(\bar{y}_{n})$ is a mono-dimensional Gaussian density $\mathcal{N}\left(\mu^2_{ij}, \sqrt{\sigma^2_{ij}}\right)$.

5.2.   Restoration of hidden process in models

In this work, we look at the MPM approach for reconstructing the hidden sequence $\textbf{X} = \textbf{x}$.

The MPM estimator for the PMC model consists of maximizing a posteriori marginal probability for each $x_n$ as follows:

where $\chi_{n}(\omega_i) = p(x_n|\textbf{y})$. Formally, this probability is calculated using the "conditional forward probabilities" $\alpha_{n}(\omega_i) = p(x_n = \omega_i|y_1, \ldots, y_n)$ and the "conditional backward probabilities" $\beta_{n}(\omega_i) = \frac{p(y_{n+1}, \ldots, y_N|x_n = \omega_i, y_n)} {p(y_{n+1}, \ldots, y_N|y_1, \ldots, y_n)}$ to avoid the numerical problems encountered by the same probabilities in the classical case. These conditional probabilities can be given using the following recursion:

● Initiation phase: For $i = 1, \ldots, K$.

(1) Forward

(2) Backward

● Induction phase: For $i = 1, \ldots, K$ and $n = 1, \ldots, N-1$.

(1) Forward

(2) Backward

with $p(x_{n+1}, y_{n+1}|x_n, y_n) = \frac{p(x_n, y_n, x_{n+1}, y_{n+1})}{p(x_n, y_n)} = \frac{p(x_n, x_{n+1})f_{x_n, x_{n+1}}(y_n, y_{n+1})}{ \sum\limits_{x_{n+1}}p(x_n, x_{n+1})f_{x_n, x_{n+1}}(y_n)}$.

It is easy to demonstrate that

Besides this, the MPM estimator under the BPMC model can be calculated as follows.

where $\chi^*_{n}(\omega_i) = p(x_n|\textbf{y}, \bar{\textbf{y}})$. This probability can be achieved using the following conditional forward and backward probabilities $\alpha^*_{n}(\omega_i) = p(x_n = \omega_i|y_1, \bar{y}_1, \ldots, y_n, \bar{y}_n)$ and $\beta^*_{n}(\omega_i) = \frac{p(y_{n+1}, \bar{y}_{n+1}, \ldots, y_N, \bar{y}_N|x_n = \omega_i, y_n, \bar{y}_n)} {p(y_{n+1}, \bar{y}_{n+1}, \ldots, y_N, \bar{y}_N|y_1, \bar{y}_1, \ldots, y_n, \bar{y}_n)}$. These conditional probabilities can be given using the following recursion:

● Initiation phase: For $i = 1, \ldots, K$.

(1) Forward

(2) Backward

● Induction phase: For $i = 1, \ldots, K$ and $n = 1, \ldots, N-1$.

(1) Forward

(2) Backward

with $p(x_{n+1}, y_{n+1}, \bar{y}_{n+1}|x_n, y_n, \bar{y}_n) = \frac{p(x_n, x_{n+1})f_{x_n, x_{n+1}}(y_n, y_{n+1})f_{x_n, x_{n+1}}(\bar{y}_{n+1})}{ \sum\limits_{x_{n+1}}p(x_n, x_{n+1})f_{x_n, x_{n+1}}(y_n)}$.

Considering the foregoing considerations, the distribution of $X_n$ conditionally on $(\textbf{Y} = \textbf{y}, \bar{\textbf{Y}} = \bar{\textbf{y}})$ is given by

This gives

5.3.   Restoration problem results

In this section, we attempt to show the importance of the BPMC model compared to the PMC model. For that, we treat an example by considering $N = 100000$ and $K = 2$ (i.e., $\Omega = \{\omega_1, \omega_2\}$). We shall investigate the models' performance by taking two factors into account, two different chain homogeneities, noted $H_1$ and $H_2$, and different noise parameters:

● Factor chain:

- Case $H_1$: $p(\omega_1, \omega_1) = p(\omega_2, \omega_2) = 0.48$ and $p(\omega_1, \omega_2) = p(\omega_2, \omega_1) = 0.02$.

- Case $H_2$: $p(\omega_1, \omega_1) = p(\omega_2, \omega_2) = 0.25$ and $p(\omega_1, \omega_2) = p(\omega_2, \omega_1) = 0.25$.

● Factor noise:

- Parameters of density $f_{\omega_1, \omega_1}$: $\mu^1_{\omega_1, \omega_1} = -3$, $\mu^2_{\omega_1, \omega_1} = -3$, $\sigma^1_{\omega_1, \omega_1} = 14$, $\sigma^2_{\omega_1, \omega_1} = 14$, and $\sigma^{12}_{\omega_1, \omega_1} = 0.1$.

- Parameters of density $f_{\omega_1, \omega_2}$: $\mu^1_{\omega_1, \omega_2} = 5$, $\mu^2_{\omega_1, \omega_2} = -5$, $\sigma^1_{\omega_1, \omega_2} = 11$, $\sigma^2_{\omega_1, \omega_2} = 9$, and $\sigma^{12}_{\omega_1, \omega_2} = 0.9$.

- Parameters of density $f_{\omega_2, \omega_1}$: $\mu^1_{\omega_2, \omega_1} = -3$, $\mu^2_{\omega_2, \omega_1} = 3$, $\sigma^1_{\omega_2, \omega_1} = 9$, $\sigma^2_{\omega_2, \omega_1} = 11$, and $\sigma^{12}_{\omega_2, \omega_1} = 0.1$.

- Parameters of density $f_{\omega_2, \omega_2}$: $\mu^1_{\omega_2, \omega_2} = 5$, $\mu^2_{\omega_2, \omega_2} = 5$, $\sigma^1_{\omega_2, \omega_2} = 18$, $\sigma^2_{\omega_2, \omega_2} = 18$, and $\sigma^{12}_{\omega_2, \omega_2} = 0.1$.

The following parameter is used to determine noise correlation influence on the error rate

$\rho = \sum\limits_{i, j = 1}^{K}\left(\frac{\rho_{\omega_i, \omega_j}}{\rho_{\omega_i^{'}, \omega_j^{'}}}\right)^2$ with $\rho_{\omega_i, \omega_j} > \rho_{\omega_i^{'}, \omega_j^{'}}$, where $\rho_{\omega_i, \omega_j} = \frac{\sigma^{12}_{\omega_i, \omega_j}}{\sqrt{\sigma^1_{\omega_i, \omega_j}}\sqrt{\sigma^2_{\omega_i, \omega_j}}}, \; \; {\mbox{for}}\; \; \; i, j = 1, \ldots, K$.

We consider four sets of $\rho$ and the covariance values associated with these four simulations are shown in Table 2.

We specify the misclassification rates computed by the MPM method in Table 3 after calculating the means from $300$ independent experiments.

According to these results, for chain in Case $H_1$, the misclassification rates are more important when the coefficient $\rho$ is small (Experiment $4$), and becomes progressively smaller when the coefficient $\rho$ decreases. It can be seen from all experiments that the BPMC model-based MPM restorations often work better than the PMC model-based ones. For the case $H_2$ of chain, the misclassification rates are consistently more important for both models. In this situation, the PMC model is preferable to the BPMC model. Other simulations have been done, showing that chain homogeneity is an important element and it can be used during the restoration phase when the data is complete.

6.   PMC and BPMC parameters estimation

The parameter estimation problem from incomplete data for the PMC is considerably more complex than for the HMM. Considering the fact that we are dealing with unsupervised segmentation, all of the parameters noted $\theta = (\theta_1, \ldots, \theta_l)$ must be estimated. We are interested, in this section, in estimating $K^2$ probabilities $p(\omega_i, \omega_j)$, $2K^2$ parameters of $\mu_{ij}$ and $3K^2$ parameters of $\Sigma_{ij}$. Considering that the log-likelihood of the models under examination cannot be maximized analytically, the methods follow the simulation according to posterior distributions used for parameter estimation. We use here ICE and SEM algorithms.

6.1.   Estimation using ICE algorithm

6.1.1.   ICE algorithm according to PMC model

The ICE algorithm under the PMC model is an iterative procedure, working as follows:

● initialize $\theta_i^0$;

● $\theta_i^{q+1}$ is calculated using $\theta_i^{q+1} = E_{\theta_i^{q}}[\hat{\theta}_i(\textbf{X}, \textbf{Y})|\textbf{Y} = \textbf{y}]$ if this expectation is computable, or $\theta_i^{q+1} = \frac{1}{T} \sum\limits_{t = 1}^T\hat{\theta}_i(\textbf{x}^i, \textbf{y})$ if the expectation above is not computable, where $\textbf{x}^1, \ldots, \textbf{x}^T$ are simulated according to $p(\textbf{x}|\textbf{y}, \theta^q)$.

Returning to our problem, we can take the following estimator for $p(\omega_i, \omega_j)$:

The conditional expectation of this estimator, at iteration $(q+1)$, can be calculated and gives

where $\psi^{(q)}_{n}(\omega_i, \omega_j) = p(x_{n} = \omega_i, x_{n+1} = \omega_j|\textbf{y})$ is a joint a posteriori probability of two consecutive classes, and can be calculated using

For the parameters $\mu_{ij}$ and $\Sigma_{ij}$, we can choose the following estimators:

However, we cannot calculate the conditional expectation of these estimators. In this situation, we employ a stochastic strategy based on the simulation of $\textbf{X}$ according to its a posteriori distribution, which is a nonstationary Markov chain, as we demonstrated previously:

with an initial distribution $\chi_1(\omega_i)$ and a transition matrix

The process for applying the ICE algorithm to the PMC model is then presented in the following manner.

● Iteration $(q = 0)$, for $i, j = 1, \ldots, K$.

(1) Initialization of the algorithm with $p^{(q)}(\omega_i, \omega_j), \mu_{ij}^{(q)}, \Sigma_{ij}^{(q)}$.

(2) Calculation of probabilities $\alpha^{(q)}_{n}(\omega_i)$ according to (5.8) and (5.10).

(3) Calculation of probabilities $\beta^{(q)}_{n}(\omega_i)$ according to (5.9) and (5.11).

(4) Deduction of probabilities $\chi^{(q)}_{n}(\omega_i)$ according to (5.12).

(5) Deduction of probabilities $\psi^{(q)}_{n}(\omega_i, \omega_j)$ according to (6.3).

(6) Simulation of $\mathbf{x}^{(q), 1}, \ldots, \mathbf{x}^{(q), T}$ according to (6.6).

● Iteration $(q+1)$, for $i, j = 1, \ldots, K$.

(1) Calculation of $p^{(q+1)}(\omega_i, \omega_j)$ according to (6.2).

(2) Calculation of $\mu^{(q+1)}_{ij}$ and $\Sigma^{(q+1)}_{ij}$ by

where $\mu^{t}_{ij}$ and $\Sigma^{t}_{ij}$ are calculated according to (6.4) and (6.5).

● Iterations are stopped if the algorithm converges; otherwise, the preceding step is repeated.

6.1.2.   ICE algorithm according to BPMC model

The ICE algorithm follows this method:

● Initialize $\theta_i^0$;

● $\theta_i^{q+1}$ is calculated using $\theta_i^{q+1} = E_{\theta_i^{q}}[\hat{\theta}_i(\textbf{X}, \textbf{Y}, \bar{\textbf{Y}})|\textbf{Y} = \textbf{y}, \bar{\textbf{Y}} = \bar{\textbf{y}}]$ if this expectation is computable, or $\theta_i^{q+1} = \frac{1}{T} \sum\limits_{t = 1}^T\hat{\theta}_i(\textbf{x}^i, \textbf{y}, \bar{\textbf{y}})$ if the expectation above is not computable, where $\textbf{x}^1, \ldots, \textbf{x}^T$ are simulated according to $p(\textbf{x}|\textbf{y}, \bar{\textbf{y}}, \theta^q)$.

Here, we take the same estimator for $p(\omega_i, \omega_j)$, as given in (6.1). The conditional expectation of this estimator, at iteration $(q+1)$, can be calculated and gives

where $\psi^{*(q)}_{n}(\omega_i, \omega_j) = p(X_{n} = \omega_i, x_{n+1} = \omega_j|\textbf{y}, \bar{\textbf{y}})$ is a joint a posteriori probability of two consecutive classes, and it can be shown that

For the parameters $\mu_{ij}$ and $\Sigma_{ij}$, we choose the same estimators as in (6.4) and (6.5).

Here, we simulate the process $\textbf{X}$ according to its a posteriori distribution under BPMC model $p(\textbf{x}|\textbf{y}, \bar{\textbf{y}})$, where

with an initial distribution $\chi^*_1(\omega_i)$ and a transition matrix

The ICE algorithm procedure for the BPMC model is then explained as follows.

● Iteration $(q = 0)$, for $i, j = 1, \ldots, K$.

(1) Initialization of the algorithm with $p^{(q)}(\omega_i, \omega_j), \mu_{ij}^{(q)}, \Sigma_{ij}^{(q)}$.

(2) Calculation of probabilities $\alpha^{*(q)}_{n}(\omega_i)$ according to (5.14) and (5.16).

(3) Calculation of probabilities $\beta^{*(q)}_{n}(\omega_i)$ according to (5.15) and (5.17).

(4) Deduction of probabilities $\chi^{*(q)}_{n}(\omega_i)$ according to (5.18).

(5) Deduction of probabilities $\psi^{*(q)}_{n}(\omega_i, \omega_j)$ according to (6.10).

(6) Simulation of $\mathbf{x}^{(q), 1}, \ldots, \mathbf{x}^{(q), T}$ according to (6.11).

● Iteration $(q+1)$, for $i, j = 1, \ldots, K$.

(1) Calculation of $p^{(q+1)}(\omega_i, \omega_j)$ according to (6.9).

(2) Calculation of $\mu^{(q+1)}_{ij}$ and $\Sigma^{(q+1)}_{ij}$ according to (6.8).

● Iterations are stopped if the algorithm converges; otherwise, the preceding step is repeated.

6.2.   Estimation using SEM algorithm

The SEM algorithm is an iterative procedure that uses stochastic drawings to estimate a sequence of model parameters from observations and realizations of $\textbf{X}$. In each iteration, we simulate the process $\textbf{X}$ according to its a posteriori distribution based on the parameters obtained in the current iteration. Then, the SEM algorithm for the models in this work proceeds as follows:

● We initialize the algorithm with $\theta_i^0$.

● At each iteration $(q)$, we simulate just one realization $\textbf{x}$ of $\textbf{X}$ according to its a posteriori distribution considering the current parameters. The parameters $\theta_i^{(q+1)}$ are then calculated.

6.2.1.   SEM algorithm according to PMC model

The SEM algorithm under the PMC model runs as follows:

● Iteration $(q = 0)$, for $i, j = 1, \ldots, K$.

(1) Initialization of the algorithm with $p^{(q)}(\omega_i, \omega_j), \mu_{ij}^{(q)}, \Sigma_{ij}^{(q)}$.

(2) Calculation of probabilities $\alpha^{(q)}_{n}(\omega_i)$ according to (5.8) and (5.10).

(3) Calculation of probabilities $\beta^{(q)}_{n}(\omega_i)$ according to (5.9) and (5.11).

(4) Deduction of probabilities $\chi^{(q)}_{n}(\omega_i)$ according to (5.12).

(5) Deduction of probabilities $\psi^{(q)}_{n}(\omega_i, \omega_j)$ according to (6.3).

(6) Simulation of $\mathbf{x}^{(q)}$ according to (6.6).

● Iteration $(q+1)$, for $i, j = 1, \ldots, K$.

(1) Calculation of $p^{(q+1)}(\omega_i, \omega_j)$ by

(2) Calculation of $\mu^{(q+1)}_{ij}$ and $\Sigma^{(q+1)}_{ij}$ by

● Iterations are stopped if the algorithm converges; otherwise, the preceding step is repeated.

6.2.2.   SEM algorithm according to BPMC model

The SEM algorithm under the BPMC model runs as follows:

● Iteration $(q = 0)$, for $i, j = 1, \ldots, K$.

● Initialization of the algorithm with $p^{(q)}(\omega_i, \omega_j), \mu_{ij}^{(q)}, \Sigma_{ij}^{(q)}$.

(1) Calculation of probabilities $\alpha^{*(q)}_{n}(\omega_i)$ according to (5.14) and (5.16).

(2) Calculation of probabilities $\beta^{*(q)}_{n}(\omega_i)$ according to (5.15) and (5.17).

(3) Deduction of probabilities $\chi^{*(q)}_{n}(\omega_i)$ according to (5.18).

(4) Deduction of probabilities $\psi^{*(q)}_{n}(\omega_i, \omega_j)$ according to (6.10).

(5) Simulation of $\mathbf{x}^{(q)}$ according to (6.11).

● Iteration $(q+1)$, for $i, j = 1, \ldots, K$.

(1) Calculation of $p^{(q+1)}(\omega_i, \omega_j)$ by

(2) Calculation of $\mu^{(q+1)}_{ij}$ and $\Sigma^{(q+1)}_{ij}$ by

● Iterations are stopped if the algorithm converges; otherwise, the preceding step is repeated.

7.   Results of images segmentation

In the present section, we show several results on the use of PMC and BPMC models for gray image segmentation in the Gaussian context. First, we present a comparative study of the results obtained from both PMC and BPMC models, concerning the misclassification rates by the Bayesian MPM algorithm, PSNR index, and the parameters estimation obtained by ICE and SEM algorithms. Second, we illustrate the segmentation results obtained on simulated images corrupted with some correlated noise. These tests give an idea about the behavior of the ICE and SEM methods for numerical and visual results.

7.1.   Numerical results

We segment some noisy simulated images with correlated noise using an unsupervised method for evaluating the endurance of the BPMC model in comparison to the PMC model, as well as the consequences when the data is not PMC or HMC-suited. For PMC model segmentation methods, the modeling of the image is done via an HPS to convert the bidimensional collection of pixels as a mono-dimensional chain, and then reorganize the image for segmentation using an inverse HPS. This transformation gives us a stochastic process with a very complex structure. This difficulty can be seen in a few papers ^[2,7,14].

In this context, the observation process representing the correlated image was obtained as follows: For each $x_s$, we simulate the variable $N_s$ according to a Gaussian distribution $\mathcal{N}(\mu_{x_s}, \sigma^2_{x_s})$, then we take

where $S$ is the set of pixels in the image in its bidimensional form, where $s_i$ is a neighbor of pixel $s$ in a neighborhood of four nearest neighbors.

We study the performance of the models presented in this work on two classes of images of different sizes, which have been noised by the previous method. However, for BPMC model segmentation methods, the sequence of pixels has been obtained via a "line by line" proceeding. We use the same form of $Y_s$ as above, and for the second element of the observable process, we set $\bar{Y}_s$ as an average of two observations for pixels neighboring to $s$ in the image but not in the chain.

The sequence of pixels collected by converting a bidimensional image to a chain of one dimension is designated by $(s_i)_{1\leq i \leq N}$. The two realizations of the stochastic processes $\mathbf{X}$ and $\mathbf{Y}$ are respectively. $\textbf{x} = (x_1, \ldots, x_n, \ldots, x_N)$ and $\textbf{y} = (y_1, \ldots, y_n, \ldots, y_N)$, where $x_i = x_{s_i}$ and $y_i = y_{s_i}$. Under these assumptions, both processes $\textbf{Z} = (\textbf{X}, \textbf{Y})$ and $(\textbf{Z}, \bar{\textbf{Y}})$ have fairly complex structures, and we can see that the distribution of $\textbf{Y}$ conditionally to $\textbf{X} = \textbf{x}$ is not necessarily a Markov; however, we have segmented the noisy image using an unsupervised method based on MPM restoration under PMC and BPMC models, where all parameters are estimated with ICE and SEM algorithms. The question is therefore whether using the BPMC model instead of the PMC model in such a context can improve segmentation results. To give some experiences, we consider two series of parameters $a$, $\mu_{x_s}$, and $\sigma_{x_s}$. Based on the parameters shown in Table 4, the original images and their noisy versions are reported in Figure 2.

We performed six segmentations, and we note that the complete data (C-D) means the original image and its noised version.

● Segmentation by MPM based on parameters of PMC model obtained from the C-D.

● Segmentation by MPM based on parameters of BPMC model obtained from the C-D.

● Segmentation by MPM based on parameters of PMC model obtained with ICE algorithm.

● Segmentation by MPM based on parameters of PMC model obtained with SEM algorithm.

● Segmentation by MPM based on parameters of BPMC model obtained with ICE algorithm.

● Segmentation by MPM based on parameters of BPMC model obtained with SEM algorithm.

We calculated some evaluated criteria to confirm the visual obtained results. The comparison of the algorithms is performed using quality measures: the error rate and the PSNR index (Table 5).

To start, we can make a number of important observations. Columns $1$ and $2$ can be seen as a kind of reference for the two models presented in this work, on the grounds that the complete data is used in the estimations. Comparing the results obtained by the estimation algorithms with those calculated by the complete data, one can notice that both the ICE and SEM algorithms implemented under the BPMC model shows favorable results indicating that misclassification rates and PSNR indexes are similar for both noise parameter series, except in the case of $\verb"Image4"$ correlated by noise with the parameter of type series $1$.

The first conclusion of interest in image segmentation is that the SEM algorithm behaves correctly in the situation under consideration. The comparison of BPMC and PMC models is the second conclusion. It is known that if the data did not fit either an HMC or a PMC, the PMC model-based MPM with HPS might have given better results than another model. This is not the case, since according to the findings and other experiments we have executed, regardless of the presence of strong enough noise, BPMC model-based MPM restorations surpass PMC model-based ones. Accordingly, using MPM under the BPMC model, ICE and SEM algorithms provide strong unsupervised segmentation methods that may substantially beat PMC model-based ones with HPS.

Also, in Table 6, parameters are also found using complete data and the parameters estimated from ICE and SEM algorithms based on BPMC and PMC models. Regarding these estimations by both algorithms, one can perceive that most of the time the estimates are similar by comparing with those obtained by C-D for PMC and BPMC models in the majority of situations of segmented images. With respect to this, the takeaway is that the BPMC model, in its current form, provides another choice that can compete with and likely outperform nearly any other model in the circumstances at hand.

7.2.   Visual results

The visual results of the unsupervised segmentation using both BPMC and PMC models are discussed here. The latter have been initialized using a two-class segmentation done by the $K$-means algorithm to generate the initial configuration of the process $X$.

The images (Figure 2) to be segmented are chosen based on the homogeneity of the objects. $\verb"Images 1"$ and $\verb"2"$ show extremely fine details, however $\verb"Images 3"$ and $\verb"4"$ are not particularly homogeneous and missing fine details. Thus, it is easy to demonstrate a significant advantage of the BPMC model over the PMC model, with HPS in the context of image segmentation. In fact, regarding Table 5, the rating is based on the average rate of each algorithm across eight experiments. Under the BPMC model, ICE was $5.25 \%$ and SEM was $2.75 \%$, while under the PMC model, ICE was $8.57 \%$ and SEM was $5.19 \%$. The PSNR index average results in the following rankings: ICE ($37.39$), SEM ($37.51$) under the BPMC model and ICE ($37.36$), and SEM ($35.84$) under the PMC model.

We can also visually observe that, depending on the segmentation phase required in these experiments, the error rate is not always the most important criterion. In fact, in $\verb"Image 4"$, the details of the zebra (Figure 3) are better restored and there are fewer spurious whites and blacks on the background in the case of the BPMC model, whereas this is not the case for the PMC model (Figure 4). However, the error rates of the two models are very close.

Nonetheless, we obtained the following conclusions:

(1) The introduction of a second component in the observation process makes the BPMC model a good rival to the PMC model with HPS in all numerical (Table 5) and visual scenarios (Figure 3) in the context of image segmentation.

(2) The SEM algorithm performs better than the ICE algorithm, especially when working on the BPMC model (Table 5).

Remark 7.1. The HPS is a space-filling curve used to read pixels during image processing. This sweep takes advantage of the two-dimensional locality of pixels. It should be noted that the introduction of such a scan is only conceivable in the case of $2^n\times2^n$ images, which is why we limit our study to the PMC model. However, the advantage of the proposed model BPMC is that it is capable of handling images of any size.

7.3.   Practical implications of proposed approach

An unsupervised approach is particularly useful when a smaller amount of or no labeled data is available for training. The proposed approach has shown superior results compared to the classical PMC model and can be potentially utilized for various applications. Particularly, the suggested model can be used for the segmentation of radar images, textured images, and medical images, because these types of images have an important correlation of noise.

8.   Conclusions and perspectives

We proposed an unsupervised approach for restoring hidden data by applying the PMC model. This study contributes to provide a 2D observed process for the latter model. The novel BPMC model is proposed first to compete with the classical PMC model for image segmentation and then to solve the problem of correlation between noise and the difficulty of two-dimensional pixel locality during one-dimensional image transformation. The parameter estimation methods described for this model are applicable to Gaussian and possibly correlated noise. We initially found that when the data follows a PMC structure, Bayesian restorations according to the BPMC model beat those based on the PMC model. Likewise, the parameter estimating methods described for the BPMC model are extremely efficient, as restorations based on complete parameters are comparable to restorations using approximated values. Next, we executed a series of tests combining the use of unsupervised segmentation of images via an HPS in the case of the PMC model and the insertion of a second element into the observed process in the case of the BPMC model. The assessment of interest in this study is based on two points. In the first step, we show that the BPMC model works best when the stochastic process from a noisy class image with correlated noise is very complex (which is neither PMC nor HMC). Several experiments in this study show that the BPMC model always gives better results compared to the PMC model, and the difference can even be quite significant in some situations. So, one way to generate complex data is to use a noisy image in a correlative version and an HPS or line by line technique. In the second step, we proved that estimation methods corresponding to the proposed model can compete with classical methods based on the PMC model with high efficiency and are therefore of interest for the image segmentation problem.

As a general conclusion, we can confirm that the proposed model, where we introduced a second component in the observed process, offers an interesting opportunity compared to the PMC model with an HPS. Basically, the composition of the second component of $\textbf{Y}$ based on the neighborhood of four nearest neighbors of each pixel provides the same role as the HPS, while on the other hand, the BPMC-based unsupervised restorations with parameter estimation methods surpass those based on the PMC model.

As perspectives for further work, we may apply the concept used here to the triplet Markovian model described in ^[40,41] in the case of Gaussian noise ^[4], in a way that the observed process will be two-dimensional and the resulting model will therefore be a bidimensional triplet Markov chain. Of course, the same concept could be used for the segmentation of fuzzy data. It consists of considering a fuzzy HMC model with bidimensional observations as presented in ^[42] and applying it to medical images. In the proposed model, both spatial correlation and locality details about pixels should be considered at the same time. Another potential direction is concerning types of noise, such as applying non-Gaussian noise to image segmentation using bidimensional Markov models.

Author contributions

A. Joumad: Conceptualization, Data curation, Writing-original draft; A. El Moutaouakkil: Conceptualization, Formal analysis, Writing-original draft; A. Nasroallah: Data curation, Formal analysis, Methodology; O. Boutkhoum: Methodology, Project administration, Software; Mejdl Safran: Funding acquisition, Investigation, Visualization; Sultan Alfarhood: Investigation, Project administration, Visualization; Imran Ashraf: Supervision, Validation, Writing-review & editing. All authors have read and approved the final version of the manuscript for publication.

Acknowledgements

The authors extend their appreciation to King Saud University for funding this research through Researchers Supporting Project Number (RSPD2024R890), King Saud University, Riyadh, Saudi Arabia.

Conflict of interest

The authors declare there is no conflict of interest.