Hyperspectral Unmixing Based on Clustered Multitask Networks

TL;DR

This paper introduces a novel hyperspectral unmixing algorithm that clusters images and employs a distributed optimization network, improving accuracy over existing methods.

Contribution

It proposes a new clustering-based, distributed optimization approach for spectral unmixing in hyperspectral images, integrating fuzzy c-means and diffusion LMS strategies.

Findings

Enhanced unmixing accuracy demonstrated in simulations

Outperforms traditional NMF-based methods

Effective clustering improves spectral signature estimation

Abstract

Hyperspectral remote sensing is a prominent research topic in data processing. Most of the spectral unmixing algorithms are developed by adopting the linear mixing models. Nonnegative matrix factorization (NMF) and its developments are used widely for estimation of signatures and fractional abundances in the SU problem. Sparsity constraints was added to NMF, and was regularized by $ L_ {q} $ norm. In this paper, at first hyperspectral images are clustered by fuzzy c- means method, and then a new algorithm based on sparsity constrained distributed optimization is used for spectral unmixing. In the proposed algorithm, a network including clusters is employed. Each pixel in the hyperspectral images considered as a node in this network. The proposed algorithm is optimized with diffusion LMS strategy, and then the update equations for fractional abundance and signature matrices are obtained. Simulation results based on defined performance metrics illustrate advantage of the proposed algorithm in spectral unmixing of hyperspectral data compared with other methods.

Figures

Tables

Table 1. Table 1 : The S ​ A ​ D 𝑆 𝐴 𝐷 SAD performance metric of five algorithms on AVIRIS Cuprite dataset, using VCA initialization.

materials	$L_{1/2}$-NMF	Dist.	S. Dist.	TV-RSNMF	Proposed
Sphene	0.2143	0.1561	0.1673	0.1583	0.1574
Nontronite	0.2518	0.1944	0.1743	0.1803	0.1711
KaolinSmect #1	0.1653	0.2370	0.1741	0.1731	0.1702
Montmorillonite	0.2318	0.3571	0.2103	0.2159	0.2248
Chalcedony	0.1995	0.1603	0.1653	0.1588	0.1437
KaolinSmect #2	0.2542	0.2873	0.2608	0.2576	0.2596
Alunite	0.3458	0.3813	0.2369	0.2551	0.2417
Buddingtonite	0.1693	0.2514	0.1953	0.2034	0.1643
Muscovite	0.1584	0.4682	0.1537	0.1563	0.1575
Andradite #1	0.3361	0.2132	0.2425	0.2392	0.2337
Dumortierite	0.2453	0.3381	0.2639	0.2686	0.2519
Andradite #2	0.3829	0.3711	0.2854	0.3136	0.2472
rmsSAD	0.2562	0.2998	0.2153	0.2207	0.2064

Equations

$${\mathbf{y}}_k=\mathbf{A}{\mathbf{s}}_k+{\mathbf{v}}_k$$

$$\sum _{n=1}^c{\mathbf{s}}_k(n)=1$$

$${\mathbf{s}}_k(n)\ge 0,n=1,...,c$$

$$J^{global}({\mathbf{s}}_1,{\mathbf{s}}_2,...,{\mathbf{s}}_N)=\sum _{k=1}^N\mathbb{E}\{|{\mathbf{y}}_k-\mathbf{A}{\mathbf{s}}_k{|}^2\}$$

$$J^{local}({\mathbf{s}}_k)=\mathbb{E}\{|{\mathbf{y}}_k-\mathbf{A}{\mathbf{s}}_k{|}^2\}+\eta \sum _{l\in {\mathcal{N}}_k\cap C(k)}{\rho }_{kl}||{\mathbf{s}}_k-{\mathbf{s}}_l|{|}^2$$

$$+\lambda ||{\mathbf{s}}_k|{|}_q$$

$$\lambda =\frac{1}{\sqrt{L}}\sum _i\frac{\sqrt{N}-||{\mathbf{y}}_i|{|}_1/||{\mathbf{y}}_i|{|}_2}{\sqrt{N-1}}$$

$${\rho }_{kj}=\frac{\theta ({\mathbf{y}}_k,{\mathbf{y}}_j)}{\sum _{l\in {\mathcal{N}}_k^{-}}\theta ({\mathbf{y}}_k,{\mathbf{y}}_l)}$$

$$\theta ({\mathbf{y}}_k,{\mathbf{y}}_j)=\frac{{\mathbf{y}}_k^T{\mathbf{y}}_j}{||{\mathbf{y}}_k||||{\mathbf{y}}_j||}$$

$$\displaystyle\mathbf{s}_{k}^{i+1}=P^{+}\big{[}\mathbf{s}_{k}^{i}+\mu\mathbf{A}^{T}(\mathbf{y}_{k}-\mathbf{A}\mathbf{s}_{k})$$

$$-\mu \eta \sum _{l\in {\mathcal{N}}_k\cap C(k)}{\rho }_{kl}({\mathbf{s}}_l^i-{\mathbf{s}}_k^i)$$

$$\displaystyle-\mu\lambda\frac{\big{(}\mathbf{s}_{k}^{i}\big{)}|\mathbf{s}_{k}^{i}|^{q-2}}{||\mathbf{s}_{k}^{i}||_{q}^{q-1}}\big{]}$$

$${\mathbf{A}}^{i+1}={\mathbf{A}}^i\ast \frac{\mathbf{Y}{\mathbf{S}}^T}{\mathbf{AS}{\mathbf{S}}^T}$$

$$||J_{new}-J_{old}||<ϵ$$

Full text

**Index Terms— ** Spectral unmixing, Hyperspectral images, Sparsity constraint, LMS strategy, Remote sensing, Distributed optimization.

1 Introduction

Remote sensing is the science of obtaining information about objects or areas from a distance, by detecting the energy that is reflected from the Earth. One of the noteworthy remote sensing techniques is hyperspectral imagery that encounter with mixed pixels challenge. Materials that are present in the scene is named endmembers. Each endmember in a pixel is weighted by its fractional abundance. The spectral unmixing (SU) technique is used to decompose a reflectance spectrum into a set of given endmember spectra. In the linear mixing model (LMM) of spectral unmixing, in which supposed that the recorded reflectance of a particular pixel, is linearly mixed by endmembers which exist in that pixel. If the number of endmembers that are present in the scene and its signatures, are unknown, the SU problem becomes a blind source separation (BSS) problem [1].

Nonnegative matrix factorization (NMF) [2], is a practical method of spectral unmixing, which decomposes the data into two nonnegative matrices. Recently, this basic method was developed by adding constraints, sparsity is one of the constraints for improving performance of NMF algorithm that is applied to the NMF cost function using $L_{q}$ regulaizers [1]. Regularization methods have been used to provide updating equations for signatures and abundances. Using $L_{1/2}$ regularization into NMF, which leads to an algorithm named $L_{1/2}$-NMF, has been proposed in [1], that enforces the sparsity of endmember abundances. Another approch is total variation regularized reweighted sparse NMF (TV-RSNMF) [3], that the total variation regularizer is embedded into the reweighted sparse NMF.

Recent researches have shown that the spatial autocorrelation of pixels gives beneficial information for the spectral unmixing, so in this paper, clustering methods and the distributed strategy has been used for utilization of spatial information. To do this, a clustering algorithm is applied to group all pixels into clusters across the whole hyperspectral image [4], before spectral unmixing, as is done in [5]. Then for using neighborhood information, a diffusion strategy such as least mean squares (LMS) is used because it has high stability over adaptive networks [6].

To solve a distributed problem, a general network which has been named clustered multitask network [6], is considered. The sparsity constrained distributed unmixing without clustering has been proposed in [7]. Here, we first cluster pixels using fuzzy c- means method (FCM), and then solve unmixing problem as a clustered multitask network using information of 8 neighboring pixels and sparsity constraint.

This paper is organized as follows. In section 2, we introduce the proposed method and optimize it. Section 3 provides simulation results and the last section gives conclusions.

2 Hyperspectral Unmixing Based on Clustered Multitask Networks

In this section, a new method that utilizes clustering of pixels and neighborhood information is proposed. First, we will express linear mixing model in subsection 2.1, and then distributed algorithm is solved for SU problem in 2.2.

2.1 Linear Mixing Model

To solve the SU problem, we focus on the linear mixing model (LMM). In this model, there exists a linear relation between the endmembers that weighted by their fractional abundances, in the scene. Mathematically, this model is described as:

[TABLE]

where $\mathbf{y}_{k}$ is an $L\times 1$ observed data vector, $\mathbf{A}$ is the $L\times c$ signature matrix, $\mathbf{s}_{k}=[s_{k}(1),s_{k}(2),..,s_{k}(c)]^{T}$ is the $c\times 1$ fractional abundance vector and $\mathbf{v}_{k}$ is assumed as a $L\times 1$ additive noise vector of $k$-th pixel of the image, when $c$, $L$ and $N$ denote the number of endmembers, bands and pixels, respectively.

In the SU problem, fractional abundance vectors have two constraints in each pixel, abundance sum to one constraint (ASC) and abundance nonnegativity constraint (ANC) [8], which are as follows, for $c$ endmembers in a scene.

[TABLE]

[TABLE]

Where $\mathbf{s}_{k}(n)$ is the fractional abundance of the $n$-th endmember in the $k$-th pixel of the image. Note that, in a BSS problem, only the observed vector is known and determination of two other matrices is our purpose.

2.2 Distributed Cost Functions and Optimization

Primarily, the fuzzy c-means clustering is adopted on the dataset [9]. Then, as explained in [6], three types of networks containing single task, multitask and clustered multitask networks are supposed. First, $N$ nodes are considered in a clustered multitask network and a optimum vector at node $k$ is estimated. A global cost function using LMS, $J^{global}(\mathbf{s}_{k}(n))$, defined as follows:

[TABLE]

where $\mathbb{E}$ is the expectation operator. Note that, the solution determined from global cost function, need to have access to information over all nodes, but the nodes can be considered to have availability only to information of its neighbors and the nodes of in the same cluster. Thus, the local cost function is used to solve this problem.

Then, the neighborhood information is used to turn the cost function to a distributed problem. In a distributed network, relationships between neighboring nodes are used to improve accuracy. In this article, we utilize the squared Euclidean distance [6], and the $L_{q}$ regularizer for sparsity constraint that is used [1]. So, the following local cost function is defined, using LMS and adding the neighborhood and sparsity constraints:

[TABLE]

where the $\mathcal{N}_{k}$ shows nodes that are in the neighborhood of node $k$, that is in the cluster $C(k)$. $\eta>0$ denotes a regularization parameter [6], that controls the effect of neighborhood term, $\lambda$ is a scalar value that weights the sparsity function [1], and the nonnegative coefficients $\rho_{kl}$ are normalized spectral similarity which is obtained from correlation of data vectors [6]:

[TABLE]

[TABLE]

where $\mathcal{N}_{k}^{-}$ include neighbors of node $k$ except itself, and $\theta$ is computed as [6]:

[TABLE]

Now, minimizing the cost function of (LABEL:eq:_10), using the iterative steepest-descent solution [10], results to:

[TABLE]

where $\mu>0$ is a step-size parameter, $i$ denotes the iteration number, and $P^{+}$ operator projects vectors onto a simplex, that adopt the ASC and ANC constraints for abundance vectors. This operator explained in [11].

Also, the spectral signatures matrix is updated similar to NMF algorithm, using multiplicative update rules [2]:

[TABLE]

A significant point in implementation of the algorithm is stopping criteria. This approach will be stopped until the maximum number of iteration ($T$), or the following stopping criteria is reached.

[TABLE]

where $J_{new}$ and $J_{old}$ are cost function values for two consecutive iterations and $\epsilon$ has been set to $10^{-8}$ in our experiments. The proposed approach is summarized in Algorithm 1.

3 Experiments and Results

In this section, for demonstration of quantitative comparison between the proposed method with different number of clusters, also between the proposed and other methods, the performance metrics such as spectral angle distance ($SAD$) and abundance angle distance ($AAD$) are used.

Primarily, the proposed algorithm has been applied on synthetic data. to generate this dataset, six signatures of USGS library have been selected randomly, using a 7$\times$7 low pass filter and containing no pure pixels. Then, the zero mean Gaussian noise with 5 different levels of SNR have been added to generated data, and performance metrics have been computed by averaging 20 Monte-Carlo runs. To choose the best number of clusters in our experiments, the SAD performance metric has been evaluated, and then according to Figure 1, the best number of clusters has been set to 6, that is equal to number of endmembers. Also, values of $\mu$ and $\eta$ has been considered equal to 0.02 and 0.1, respectively [6], and $q=1$, to gain the best results. Then the proposed algorithm and some other algorithms such as VCA-FCLS [12], NMF [2], $L_{1/2}$-NMF [1], distributed unmixing [6], sparsity constrained distributed unmixing [7] and TV-RSNMF [3], that is similar to the proposed algorithm without clustering step, has been applied on the generated synthetic dataset. The comparison of performance metrics of this 7 different methods has been shown in Figure 2, where the metrics of proposed algorithm is star-dashed line and excels other methods. Afterwards, the proposed algorithm has been applied on AVIRIS Cuprite real dataset [13]. After clustering into 12 clusters, the simulation results of spectral signatures and fractional abundances have been shown in Figure 3 and  4. Also, $SAD$ performance metric of $L_{1/2}$-NMF, distributed unmixing, sparsity constrained distributed unmixing and proposed method on the real dataset have been compared in Table 1, the results of proposed algorithm are available in the last column and has the best $rmsSAD$ value.

4 Conclusion

Spectral unmixing (SU) is a technique to characterize mixed pixels in hyperspectral images measured by remote sensors. Decomposition of pixels in the scene into their constituent materials is the goal of spectral unmixing. This paper followed two steps, including the FCM clustering of hyperspectral images and then the sparsity constrained distributed unmixing method. This new algorithm considered sparsity, clustering and neighborhood information. Simulation results on synthetic and real datasets illustrated preference of proposed approach in comparison against NMF, $L_{1/2}$-NMF, VCA-FCLS, distributed unmixing and sparsity constrained distributed unmixing methods.