{'fraction_non_alphanumeric': 0.026485204427377458, 'fraction_numerical': 0.014230856110232663, 'mean_word_length': 4.8522802379378716, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Using Spatial Correlation in Semi-Supervised Hyperspectral Unmixing under Polynomial Post-nonlinear Mixing Model

This paper presents a semi-supervised hyperspectral unmixing solution that integrate the spatial information in the abundance estimation procedure. The proposed method is applied on a nonlinear model based on polynomial postnonlinear mixing model where characterizes each pixel reflections composed of nonlinear function of pure spectral signatures added by noise. We partitioned the image to classes where contains similar materials so share the same abundance vector. The spatial correlation between pixels belonging to each class is modelled by Markov Random Field. A Bayesian framework is proposed to estimate the classes and corresponding abundance vectors alternatively. We proposed sparse Dirichlet prior for abundance vector that made it possible to use this algorithm in semi-supervised scenario where the exact involved materials are unknown. In this approach, we just need to have a large library of pure spectral signatures including the desired materials. An MCMC algorithm is used to estimate the abundance vector based on generated samples. The result of implementation on simulated data shows the prominence of proposed approach.


Introduction
Hyperspectral images are mentioned as images have been taken in a hundred of spectral bands in remote sensing area [Keshava 2002]. They have been used in various applications including agricultural and environmental monitoring [Somers 2009], mineral exploration [Settle 1993], military surveillance [Chang 2000], and so on [Liu 2014]. Due to physical limitations of imaging devices, each hyperspectral pixel is mixture of 4 one of these publicly available libraries which has been taken over 22 years covers more than 1,300 spectral signature of so many materials [Clark 2007]. Since, in each hyperspectral pixel only a small number of endmembers against the extremely large library are contributed, the abundance vector could be sparse. In this paper a sparse Dirichlet prior is proposed for abundance vector. Thus, unmixing and endmember selection from a large dictionary are executed simultaneously. We also make use of MRF to profit from spatial correlation in powerful unmixing algorithm Bayesian PPNMM and thus major improvement in unmixing accuracy is achieved. Using MRF during unmixing procedure, made it possible to classify the hyperspectral image to known number of classes wherein the abundance vector and the nonlinearity term is the same in each class.
Thus our algorithm performs not only endmember selection during unmixing process in a flexible generalization of LMM, PPNMM, but also make classification of the hyperspectral pixels.

Problem Formulation
It is important to note that due to contribution of only a small number of endmembers of an extremely large library in each hyperspectral pixel, the abundance vector could be sparse. Accordingly, in this paper we consider this case for estimating abundance vectors in a Bayesian sense. In this way, unmixing and endmember selection from a large library are performed simultaneously.
To elaborate, in a nonlinear mixing model, a hyperspectral pixel is defined as a nonlinear function of a linear mixture of endmember signatures affected by noise term as where is an -dimensional hyperspectral pixel, is the spectral signature of the th endmember in the library , is the corresponding abundance, is the number of endmembers in the library, is a nonlinear transformation, and is additive white Gaussian noise with zero-mean and variance 2 : where denotes an × identity matrix.
Here, a second-order polynomial is employed for the nonlinear function . It has been shown that second order polynomial is an appropriate approximation for nonlinear mixing models [Nascimento 2009], since higher order terms are negligible and could be merged in noise term. A PPNMM [Altmann 2012], is given by where ʘ is the Hadamard product. This model contains both bilinear and linear models.
If = 0, the model is simplified to LMM and thus the result would be at least as good as linear ones. According to [Altmann 2012], using as a single amplitude parameter for the nonlinear term, lower complexity is achieved.
Abundances has some consideration in practise. There are two constraints known as non-negativity and sum-to-one that mean each element of abundance vector should not be smaller than zero and sum of the abundance fractions should be equal to one:

Spatial Correlation Formulation
Previous efforts exhibited that MRFs are interesting tools for modelling spatial correlation used in hyperspectral image classification and segmentation in Bayesian framework [Chen 2017]. We utilize Potts-Markov random field [Eches 2011] where ( ) is neighborhood of pixel , and \ = { , ≠ }. Using Potts-Markov model, based on Hammersly-Clifford theorem which relates the MRFs to Gibbs distribution, the pdf of random field is written as follows: where > 0 is known as granularity coefficient, ( ) is normalizing constant or  In this paper we consider 4-pixel first order neighborhood structure, so according to [Eches 2011], a fixed moderate value of is enough to avoid from trapping in a local optimum. For larger neighborhood structure, smaller value of is appropriate.

Hierarchical Bayesian Framework
In this section the likelihood function of observed pixel is computed based on PPNMM for the hyperspectral unmixing and the priors are considered to calculate the posteriors of unknown parameters. We utilize the hierarchical Bayesian solution. Our motivation is to select a proper prior for the abundance vector which leads to not using any EEA. Also we benefit the spatial correlation in nonlinear mixing model to enhance the unmixing performance.
Assuming Gaussian noise, the likelihood function of the mixed pixel is normally distributed denoted as: where is the observed hyperspectral pixel and is its corresponding class label.
Unknown parameters , , and 2 should be estimated. Due to independency between noise elements for pixels = 1, … , we have:

Prior Selection
In this section the prior distributions that have been chosen for unknown parameters = { , , , 2 } are described. In addition, a hyperparameter could be introduced and a proper prior would be assigned to it.

Pixel Class Prior:
As explained earlier, the prior pdf of the pixel class vector is a Potts-Markov random field with first-order neighborhood as [Eches 2011]. The parameter for the distribution is set to 1.1.

Abundance Prior:
Due to the two physical constraints, sum-to-one and non-negativity of the abundance vector, we propose the sparse symmetric Dirichlet distribution [Ng 2011] as a prior for the abundance vector of class as . This distribution presents a sparse behavior for < 1 and corresponds to a uniform distribution over the standard ( − 1)-simplex for = 1. The latter case is commonly used in unmixing problems [Altmann 2012]. Note that this case is limited to supervised unmixing applications in which the exact endmembers must be known. In more realistic scenarios, however, exact mixing endmembers are unknown and only a large spectral library is given and the concentration parameter plays an essential role. Here, accordingly we consider < 1 and show that this case is suitable for a wide range of applications.

Noise Variance Prior
We use a non-informative Jeffreys' prior for noise variance 2 [Altmann 2012] as follows: where ℝ + (. ) is the indicator function defined on the positive real values:

Nonlinearity Coefficient Prior
For the unknown parameter , the following priors is assigned: where 2 is a hyperparameter for which the Inverse-Gamma prior is where ( , ) are set to (1, 0.01) according to [Altmann 2012].

Derivation of Posterior Distribution
Using the likelihood of mixed pixel and assumed priors, the joint posterior distribution of all unknown variables = { , 2 } would be extracted by hierarchical structure as: which means ( | ) is proportional to the product of the likelihood function by the priors and hyperpriors.
As seen in (15) [Robert 2004] to these samples.

Metropolis-Within-Gibbs Sampler
In this section we investigate the Metropolis-within-Gibbs sampler to generate samples iteratively according to distributions ( | ) wherein there exist unknown parameters { , , , 2 } and one unknown hyper-parameter 2 .

Conditional distribution of
For each pixel , the Bayes theorem yields As seen this distribution is too complex to be directly sampled. Here, we use the MCMC sampler to generate samples distributed according to (17). As mentioned earlier, we choose Metropolis-within-Gibbs sampler in which new samples are generated by Gaussian random walk procedure [Roberts 1996].

Conditional distribution of
Using Bayes theorem, the conditional distribution of class label is expressed by the probabilities The derived posterior distribution of the class labels of all pixels in (18)

Conditional distribution of
Using (7) and the Gaussian prior for nonlinear coefficient , the conditional distribution is obtained Gaussian too, as follows  (21) where ( ) = ( ) ʘ ( ), which ʘ is Hadamard product, and is number of pixels belong to class label . Consequently, sampling from (19) is done easily.

Conditional distribution of 2
The full-conditional pdf of noise variance is expressed as following based on Bayes theorem: Due to choosing Jeffrey's prior, the conditional distribution is described as which it is a known distribution that can easily be sampled.

Conditional distribution of 2
Considering (15), straightforward computations leads one to the following distribution from which it is easy to sample.
To end, using computed conditional distribution of unknown parameters and hyperparameter, a summary of algorithm can be provided as Alg. 2.
Algorithm 2 Gibbs sampler proposed for hyperspectral unmixing using spatial correlation in PPNMM  (19) 11. Sample 2( ) according to (23) 12. Sample 2( ) according to (24) Following the explained procedure in Alg. 2, the samples for all unknown parameters would be generated for which an MMSE estimator could be used to compute the sample mean of them as the estimates.

Experimental Results
To evaluate the performance of proposed unmixing algorithm, our experiments were performed on both synthetic and real hyperspectral image.

Simulated Data
In first evaluation strategy, three 25 × 25 pixel hyperspectral images with = 3 classes were generated. These images 1 , 2 and 3 were generated by different mixing models including LMM, GBM and PPNMM to evaluate the robustness of algorithms accuracy to different mixing models. The nonlinear parameter = 0.1 was used for PPNMM synthesized image, and similarly the nonlinear parameter = [ 0.5 0.1 0.3 ] was established for GBM synthesized image. For each of = 3 label classes, a unique abundance vector was used according to Table 1. An additive Gaussian noise with variance 2 = 0.001 corrupted the synthetic images. This variance corresponds to a reasonable signal to noise ratio around 15 dB in our problem. We select 8 endmembers randomly from the USGS library [Clark 2007] and make our own library. We selected three endmembers to contribute making mixed pixel while five other endmembers were just in library to simulate a semi-supervised scenario. Note that in this case we are not aware of neither the number of endmembers nor the associated ones in the mixing process.
The spectral signature of these materials are illustrated in Fig. 2.

Figure 2 Spectral signature of three endmembers from USGS [Clark 2007]
A random label map which generated by MRF with = 1.1 was assigned to each images. We also choose = 0.2 making a sparse distribution. Simulations run with 5000 MCMC and 500 burn-in iterations.  The real and estimated class labels are illustrated in Fig. 1 (a), 2 and 3 for LMM, GBM and PPNMM synthesized images respectively. The initial random class labels and estimated class map for three first iterations are shown in Fig. 1 (b) for linear mixing model. As seen, the class labels converge to real label map in initial iterations and finally goes to true label map. The results demonstrate that the algorithms estimate the pixel classes perfectly. The classification is done by unmixing simultaneously.
(26) Table 2 shows the RMSE of abundance estimates on three synthesized images for three most common unmixing algorithms. The implemented unmixing algorithms were classified to three categories: the LMM-based solutions, the GBM-based solutions and the PPNMM-based solutions. Since our aim was to have a comprehensive comparison, for each type of mixing model a corresponding appropriate unmixing algorithm was used. It's important to note that FCLS, MRF based LMM algorithm, Sub-gradient GBM, Sub-gradient PPNMM and Supervised Bayesian PPNMM algorithm were performed on a supervised manner; which means the exact three endmembers were fed to unmixing procedure. While, in our proposed method a library contains eight endmembers was provided to unmixing process. We evaluate the Bayesian PPNMM algorithm in semi-supervised scenario to only compare the result with our algorithm. The results show that if the materials in library get increased compared to true endmembers, the supervised algorithms would dropped significantly in estimation accuracy. However, our method is not so sensitive to the size of library. As seen in Table 2, best results are acquired by proposed algorithm. Comparison between MRF based unmixing algorithms and the others show that using spatial correlation enhance the unmixing accuracy significantly. Also, using MRF in PPNMM results smaller estimation error compared to the LMM one even in 1 where the synthesized image is made under linear mixing model. In Table 3 the REs of reconstructed image for three synthesized images were presented. We computed the RE respect to the ideal synthesized image, not to the noisy one. Because if the algorithm had been sensitive to noise, the RE for noisy image would be smaller than the ideal image. As seen, the RE for the MRF-Bayesian is around ten times smaller than other algorithms.

Real Data
The second part of evaluation was performed on real data collected by the Airborne We select 93 endmembers of the USGS mineral library as our own library, whereas for supervised MRF-LMM algorithm the endmembers are extracted from VCA algorithm [Nascimento 2005] with = 14 as pointed in [Eches 2011]. The algorithm was run with 5000 MCMC iterations and 500 burn-in iterations. We set the number of classes equal to = 5 according to [Eches 2011]. The resulted estimated abundance vectors are illustrated in Fig. 6. As seen, most of the abundance values for each class are near zero, so the sparsity constraint is valid.
The abundance values for the most significant endmembers are also shown in Table 4 for 5 classes. It's interesting to note that if the number of classes set larger than the actual number of classes, the estimated abundance value corresponding to the same classes would be similar. This result is understanding from Table 4 for the second and fifth classes.  We also showed the abundance map for six most significant endmembers in Fig. 8 and the classification map was plotted in Fig. 9. Figure 8 shows that the proportion of materials in disjoint regions are completely different, which this information helps to find the classes of all pixels. Looking at the label maps in Fig. 9, it is obvious that although  The reconstruction error is computed for two best algorithms of Table 2 and reported in Table 5. An interesting result was achieved; both algorithms had the same accuracy in the sense of reconstruction error. However there is no measure to compare the abundance estimates because for real images, the exact proportion of materials is not available. It's important to note that although both algorithms result equal error, but the evaluation condition is different; LMM based algorithm is applied in supervised manner, while the proposed algorithm is in semi-supervised scenario. The former one had the 24 exact 14 endmembers whereas the latter should select 14 endmembers from a library of 93 members. Under this rigid situation the proposed algorithm performed excellently. We repeat the algorithm for a larger hyperspectral image from Cuprite image.
The region of interest of size 200 × 200 is shown in Fig. 10. We set the number of classes = 14 according to [Eches 2011]. We use the reduced mineral USGS library with 93 endmembers. The estimated label map is shown in Fig. 11. As seen the algorithm could labeled the homogenous parts of image accurately.

Conclusion
We derived a hierarchical Bayesian algorithm for unmixing of hyperspectral images based on the Polynomial Post-nonlinear Mixing Model. A Dirichlet prior was proposed for modeling the abundance vector sparsity. We set the concentration parameter in such way that the abundance pdf leads to a sparse distribution. In this way, if a huge library is given, the unmixing procedure could be done precisely and any third-party EE algorithm would not be necessary. We also perform Markov Random Fields to the nonlinear mixing model to benefit from spatial correlation to enhance the unmixing process. So hyperspectral classification could be done simultaneous to nonlinear unmixing. Due to complexity of derived posterior distribution, the MCMC method was used to estimate the posterior. Using proposed algorithm, the unmixing accuracy improved compared to LMM, GBM and simple PPNMM algorithms in the sense of RMSE and RE. Also, although we utilize the proposed algorithm in semi-supervised manner, but the nonlinear based MRF results more accurate label map in contrast to linear one.