Hyperspectral Endmember Extraction Techniques

Abstract

Hyperspectral data processing and analysis mainly plays a vital role in detection, identification, discrimination and estimation of earth surface materials. It involves atmospheric correction, dimensionality reduction, endmember extraction, spectral unmixing and classification phases. One of the ultimate aims of hyperspectral data processing and analysis is to achieve high classification accuracy. The classification accuracy of hyperspectral data most probably depends upon image-derived endmembers. Ideally, an endmember is defined as a spectrally unique, idealized and pure signature of a surface material. Extraction of consistent and desired endmember is one of the important criteria to achieve the high accuracy of hyperspectral data classification and spectral unmixing. Several methods, strategies and algorithms are proposed by various researchers to extract the endmembers from hyperspectral imagery. Most of these techniques and algorithms are significantly dependent on user-defined input parameters, and this issue is subjective because there is no standard specificity about these input parameters. This leads to inconsistencies in overall endmember extraction. To resolve the aforementioned problems, systematic, generic, robust and automated mechanism of endmember extraction is required. This chapter gives and highlights the generic approach of endmember extraction with popular algorithm limitations and challenges.

1. Introduction

The hyperspectral imaging (HSI) is a renowned technology that uses the chemical composition based spectroscopic properties of earth surface materials for their detailed analysis and exploration [1]. HSI is spectrally overdetermined and has the ability to observe the surface materials continuously across the wide range of electromagnetic spectrum, generally covering 0.4–2.5 μm, where different regions of wavelengths gives information about different material contents [2]. HSI, rather than comprehensively depending only on the spatial variations within the image (as in multispectral imaging), has moved to take the advantage of spectral variations for detailed material analysis [3]. The narrow and continuous spectral measurement provides the enriched amount of information significantly usable for target detection, material identification, material mapping, surface material classification, abundance estimation and mapping details of surface properties [2, 4]. These capabilities of HSI led to extend its areas of application and actively appealing researchers coming from varying subject domains to solve their respective problems. The popular areas, where hyperspectral imaging is providing its significance includes agriculture, forestry, hydrology, oceanography, soil analysis, land use land cover mapping and geology, etc.

Along with the ample amount of information, the hyperspectral imaging is characterized by several data processing and analysis challenges including adverse atmospheric effects, curse of dimensionality, unavailability of the ground-truth information and spectral mixing [1, 2, 3]. To extract the maximum information from this high volume HSI data with greater precision, these challenges needs to be handled efficiently. Of course, there several popular digital image processing and pattern recognition techniques are well available in the literature, but those techniques cannot be used as it is to process HSI data due to its different challenges. To process and analyses the hyperspectral data proficiently with better accuracy, there is a need to follow the generalized and application independent data processing framework ( Figure 1 ) [1].

Figure 1.

Generalized hyperspectral data processing framework.

The ultimate intention of hyperspectral image analysis is to achieve the better classification or mapping accuracy, but it is adversely affected due the two important issues. First, the unavailability of the ground truth information or lower number of training samples leads to the Hughes effect, which decreases the classification accuracy as the number of features increases. Second, the spectral mixing problem raised either due to the limited spatial resolution of the sensor or homogenously combined mixture of surface materials. This leads to the linearly combined spectral reflectance of the component materials for every mixed pixel and keep that particular pixel away from any of the class in classification. To resolve these significant obstacles of HSI data classification and mapping, the term, image derived “endmember” comes into picture and contribute with remarkable role in increasing HSI data classification accuracy [1].

2. Hyperspectral endmembers

Due to significantly enhanced spatial and spectral resolutions of the hyperspectral imaging sensors, hyperspectral endmember extraction (EE) is found to be a significant step in hyperspectral data analysis. As described in [5], “a hyperspectral endmember (also known as ‘pure pixel’) is an idealized, pure signature of a spectral class.” The pure spectral signature signifies the complete reflectance of pixel exclusively occupied by a single surface material. In case of multispectral imagery, finding endmembers within the scene is significantly difficult because of the broad spectral bands with lower spatial and spectral resolutions. Therefore the EE is not taken as an important step in multispectral data analysis and its importance has remained unnoticed. Thereafter, with the substantial improvements in hyperspectral sensors several subtle surface material contents that cannot be explored by multispectral imaging can now be explored using the hyperspectral imagery. These material contents are usually not known previously and can be only analyzed with the higher spectral resolution data. The hyperspectral endmembers are supposed to be one of these material contents. In HSI data, finding endmembers is a critical task due to their population within the data is significantly low and generally beyond the human visual perception [6]. Once the endmembers are extracted precisely from the HSI image itself, they can be used as a reference data or training data to enhance the classification accuracy or efficient spectral unmixing of the hyperspectral images [1]. These image derived endmembers are more effective as compared to the standard laboratory or field spectra, as they are not certainly recorded under the similar conditions as the satellite or airborne HSI data.

There are several EE techniques and methods are available in the literature, each having its own advantages and challenges. Depending on the EE approaches, these techniques are broadly classified into two categories ( Figure 2 ). First, the convex geometry based EE. Further, it is sub-classified into the Orthogonal Projection (OP) based and Simplex Volume (SV) based EE approaches. The OP based EE approaches, makes the orthogonal projection of all data samples onto a set of selected vectors and considers the data samples producing extreme (either minimal or maximal) projections with these selected vectors as a final set of endmembers. The popular OP based EE algorithms are Pixel Purity Index (PPI) [7], VCA [8] and Sequential Maximum Angle Convex Cone (SMACC) [9]. The SV based EE approaches, assumes that the simplex formed by a set of pure signatures as vertices should produce the maximum volume among all simplexes formed by the same number of signatures as vertices. The SV based EE algorithms are NFINDR [10], Simplex Growing Algorithm (SGA) [11] and Convex Cone Analysis (CCA) [12]. Second, the statistics based EE. Further, it is sub-classified into second order statistics based EE, higher order statistics based EE and Swarm Intelligence based EE. The Iterative Error Analysis (IEA) [13] and EESCA [14] techniques comes under the category of second order statistics based EE. The skewness and Independent Component Analysis EE algorithm (ICA-EEA) [6] comes under the category of higher order statistics based EE. Apart from these traditional techniques, few significant Swarm Intelligence (SI) based approaches are also discussed. The Ant Colony Optimization (ACO) [15], Particle Swarm Optimization (PSO) [16] and Artificial Bee Colony (ABC) [17] and others are the typical examples of SI based EE algorithms. SI is the new area in artificial intelligence and deals with the natural and artificial systems composed of multiple entities that are articulated using self-organizations and decentralized control. In particular, the discipline emphases on the mutual behaviors that resulted from the interactions of entities with each other & their environment. SI algorithms have capability to solve the problem of combinatorial optimization.

Figure 2.

EE techniques.

In the further portion, this chapter present and discuss these EE algorithms in step-wise manner along with few practical experiments. The Samson hyperspectral image ( Figure 3 ) is used for practical demonstrations. This image contains 156 spectral channels (having 3.13 nm spectral resolution) covering from 401 to 889 nm range of electromagnetic spectrum. Spectrally, there are three classes are available into this image, i.e., soil, tree and water [18].

Figure 3.

Samson hyperspectral image with ground truth (GT). (a) Samson. (b) GT: abundances. (c) GT: endmembers.

2.1 Pixel purity index (PPI)

The PPI algorithm [7] is a mostly used endmember extraction algorithm due to its availability into the ENVI commercial software package. Its step-by-step implementation of is never been revealed into the literature due to its propriety rights. Considering its basic idea, several endmember extraction algorithms [19, 20, 21] tried to mimic the implementation of original PPI [7] and produced the useful results. The PPI is an orthogonal projection based endmember extraction mechanism and works through stages briefly explained below:

2.1.1 Initialization

The algorithm initiates with the randomly generate “k” unit vectors, also called as “skewers,” ${\left\{{\mathit{\text{skewer}}}_k\right\}}_{k=1}^k$ where k is a supposed to be user defined large positive integer.

2.1.2 PPI count calculation

In PPI count calculation, every data sample is orthogonally projected onto all the skewers ${\left\{{\mathit{\text{skewer}}}_k\right\}}_{k=1}^k$ . Further, for every data sample r, identify those skewers onto which data sample r produces extreme (either the maximal or minimal) projection at their end points. These skewers producing extreme projections, form a set (i.e., $S_{\mathit{\text{extrema}}}\left(r\right)$ ) for that particular data sample r. Then, the PPI count of that particular data sample r is equal to the cardinality (i.e., total number of skewers within the set) of the set $S_{\mathit{\text{extrema}}}\left(r\right)$ , that is,

$$n_{\mathit{PPI}}\left(r\right)=\mid S_{\mathit{\text{extrema}}}\left(r\right)\mid$$ E1

where |A| is defined as the cardinality of set A.

2.1.3 Candidate selection

In candidate selection, the user defined value is set for threshold t to threshold the PPI count $n_{\mathit{PPI}}\left(r\right)$ for all data sample vectors to select the possible data samples for endmembers.

2.1.4 Endmember extraction

At the end, all the data samples having $n_{\mathit{PPI}}\left(r\right)\ge t$ as considered as final set of endmembers.

Figure 4 enlightens the endmembers extracted using PPI from Samson data. In this experiment, PPI is executed with 1000 iterations and 50 as a threshold value, where it has extracted seven signatures as a final set of endmembers. Out of these seven endmembers, one represents the “soil”, four represents the “trees” and two represents the “water.” One endmember signature (along with their respective spatial coordinates) for each class in the Samson image is given in Figure 4 along with ground truth spectra’s.

Figure 4.

PPI extracted endmember from Samson data.

Although PPI is mostly used endmember extraction algorithm, but initially it was not considered as a solution for endmember extraction, but it was considered as a guide for endmember extraction. The PPI is very sensitive to the value of randomly generated skewers and threshold. The need of large number of randomly generated skewers leads to the higher computations complexity and yields varying outcomes (due to randomly generated skewers) during several runs even on the same hyperspectral image. Further the undefined criteria for selecting the appropriate value of threshold put forward the need of trained user for selection the final set of endmembers [19]. To deal with these challenges, there are few attempts available into the literature that tried to improve the performance of PPI with some modifications at algorithm level. In [19, 20], the attempt is made use the concept of Virtual Dimensionality (VD) to automatically identify the number of endmembers to be identified from the image and opted the algorithm based endmember initialization to produce the appropriate initial endmember set, which minimizes the large number of runs required for PPI. In [21], multi-dimensional PPI (MDPPI) is presented as a fast alternative for PPI, which iteratively identifies the convex hull indices in lower-dimensional random projections.

2.2 Vertex component analysis (VCA)

The vertex component analysis (VCA) algorithm [8, 22, 23] is an unsupervised endmember extraction algorithm and works with the assumption that, in linear spectral mixing, every pixel signature is composed with the linear combinations of endmember spectra available within the scene. The VCA algorithm explores two facts: One, the endmembers are found to be the vertices of the simplexes and two, the affine transformation of every simplex is too simplex. The VCA algorithm initiates with the assumption that, there is a presence of endmembers within the data and iteratively projects data sample vectors on to the direction orthogonal to the subspace covered by the hyperspectral endmembers which are already identified. The new endmembers relates to the extreme (either minimal of maximal) projections. The algorithm continues to iterate till the number of endmembers are exhausted. The step-wise VCA algorithm is briefed below:

1. Suppose $r$ be the $L\times 1$ vector, where L represents the total number of spectral channels, and $m_i$ is the spectral signature of $i^{\mathit{th}}$ endmember, therefore

   $$r=\mathit{e\alpha },$$ E2

   where, $e=\left[e_1,e_2,\cdots ,e_p\right]$ and $\alpha ={\left[{\alpha }_1,{\alpha }_2,\cdots ,{\alpha }_p\right]}^T$ is the abundance portion of every endmember, and $p$ is the total number of endmembers available within the image.

2. Due to the non-negativity constraint, the endmember abundance fractions satisfy the conditions given in Eqs. (3) and (4).

   $$0\le {\alpha }_k\le 1$$ E3
   $$\sum _{k=1}^p{\alpha }_k=1$$ E4

3. Every data sample can be considered as a sample vector within the L-dimensional Euclidean space, where every channel is allocated to the one axis of space and all are mutually orthogonal. Because of the constraints in Eqs. (3) and (4), the observed sample vector $r$ belongs to the simplex having endmembers in the vertices. The detailed implementation and evaluation of VCA is given in [8].

Figure 5 gives the endmembers extracted using VCA from Samson data. In this experiment VCA is executed to extract three endmembers using nine iterations. It has extracted three endmembers with spatial coordinates (69, 28), (4, 84) and (14, 54), out of which, first signature represents the “soil,” second signature represents the “tree” and third signature represents the unknown or may be mixed pixel.

Figure 5.

VCA extracted endmember from Samson data.

The VCA algorithm is considered as a sequential implementation of the PPI algorithm. It is having lower computations cost as compared to the PPI and NFINDR EE algorithms. The usage of VCA is requires the precise knowledge about the number of endmembers to be identified. Its random initialization nature leads to the inconsistent outcomes during several executions even of the same image.

2.3 Sequential maximum angle convex cone (SMACC)

The sequential maximum angle convex cone (SMACC) [9] endmember extraction algorithm make use of a convex cone model along with non-negativity or sum-to-unity constraints to identify hyperspectral endmembers. In SMACC, the extreme projection points are used to determine a convex cone, which describes the first endmember. Further the constrained oblique projection is applied to the existing cone to extract another endmember. The convex cone is enlarged to comprise the new endmember. This procedure is revised till a projection extracts the endmember that is already available in the convex cone or till the defined number of endmembers is extracted.

In other words, SMACC algorithm first identifies the brightest pixel from the image; thereafter, it identifies the pixel mostly diverse from the brightest. Thereafter, it identifies the pixel that is most diverse as compared to the first two pixels. This procedure is revised until SMACC identifies a pixel which is already accounted for within the group of the formerly identified pixels, or until it identifies a previously defined number of endmembers.

As the convex approaches that are dependent on a simplex analysis, the number of hyperspectral endmembers is not limited by the total number of spectral bands. Although the endmembers identified by SMACC are exclusive, a one-to-one correspondence is not available among the total number of materials within the image and the total number of endmembers. SMACC extracts endmembers from pixels within the image. Every image pixel may be occupied by only single material or may be by high proportion of a single material with unique combinations of other surface materials. Every single material identified from the image is defined by a subset covering its spectral inconsistency. SMACC algorithm provides the endmember basis that describes each of these material subsets. Apart from this, SMACC also offers abundance images to identify the fractions of the total spectrally integrated pixels contributed by each constituent endmembers. Mathematically, SMACC algorithm uses the following convex cone expansion for every endmember spectra, defined as p:

$$p\left(c,i\right){=}_k^{N}R\left(c,k\right)A\left(k,j\right)$$ E5

where, $i$ indicates the pixel index, $j$ and $k$ represents the endmember indices from 1 to the expansion length N, $R$ is the matrix containing endmember signatures as columns, $c$ gives the spectral channel index and $A$ gives the matrix that comprises the proportional contribution of every individual endmember $j$ in every endmember $k$ for every pixel.

The SMACC EE tool available into the ENVI software used for EE from Samson data. The SMACC is executed to extract three endmembers. As compared to the GT in Figure 6 , SMACC has successfully extracted two endmembers signifying the “soil” and “tree” classes. The third signature belongs to the unknown class, which is incorrectly reported as an endmember by SMACC.

Figure 6.

SMACC extracted endmembers from Samson data.

The overall observation is that, SMACC provides a faster and more automated method for finding spectral endmembers, but it is more approximate and yields less precision.

2.4 NFINDR

The NFINDR [10] algorithm is basically based on the geometric properties of the convex sets. The fundamental idea of NFINDR is to identify data samples that can produce a simplex having maximum volume and these samples are considered as final endmembers. Wherein, the algorithm assumes that, in L-spectral dimensions, the L-volume contained by a simplex produced of the purest data sample is bigger as compared to any other volume produced by any other combinations of data samples. The algorithm proceeds by “inflating” a simplex within the data, initializing with the randomly selected set of pixels. For every data sample and for every endmember, the endmember is substituted with the spectra of the data sample and the volume is re-calculated. If the volume founds increased, the spectra of the new data sample substitutes that endmember. This procedure is revised till no any replacement takes place. The step-wise implementation along with the mathematical perspective is briefly given below:

Usually, the spectrum of a given data sample is supposed to be linear combinations of the endmember signatures.

$$p_{\mathit{ij}}=\sum _k{e}_{\mathit{ik}}{c}_{\mathit{kj}}+\epsilon$$ E6

where $p_{\mathit{ij}}$ indicates the ith spectral channel of jth data sample, $e_{\mathit{ik}}$ indicates the ith spectral channel of kth pure pixel, $c_{\mathit{kj}}$ denotes the proportions of spectral mixing for the jth data sample from the kth pure pixel and $\epsilon$ indicates the Gaussian error assumed to be very small. As the sample compositions in spectral mixing are expected to be in percentage, all the mixing proportions must sum to one:

$$\sum _k{c}_{\mathit{kj}}=1$$ E7

If, $c_{\mathit{kj}}\approx 1$ for any of the endmember contribution within the data sample, the other pure pixel contribution is considered nearly equal to zero and the data sample can be considered as an endmember.

2.4.1 Preprocessing

In preprocessing, the orthogonal subspace projection (OSP) and maximum noise fraction (MNF) transformations are applied to the dataset to minimize the data dimensionality one less than the total number of endmembers to be identified. It helps to determine the data sample volume in a simpler manner.

2.4.2 Volume determination

Consider E as an endmember matrix augmented by row of ones:

$$E=\left[\begin{array}{c}1\\ {\overrightarrow{e}}_1\end{array}\begin{array}{c}1\\ {\overrightarrow{e}}_2\end{array}\begin{array}{c}\cdots \\ \cdots \end{array}\begin{array}{c}1\\ {\overrightarrow{e}}_i\end{array}\right]$$ E8

where ${\overrightarrow{e}}_i$ is the column vector comprising the spectrum of endmember i.

Now, the simplex volume (V) produced with the endmember estimates is directly proportional to determinant E.

$$V\left(E\right)=\frac{1}{\left(l-1\right)!}\mathit{abs}\left(\left|E\right|\right)$$ E9

Here, $\left(l-1\right)$ represents the number of dimensions occupied by data.

2.4.3 Endmember selection

The algorithm is stopped when all the pixels are tested and no replacement takes place.

Figure 7 gives the endmembers extracted using NFINDR from Samson data. In this experiment NFINDR is executed to extract three endmembers using nine iterations. It has successfully extracted three endmembers with spatial coordinates ((69, 29), (4, 84) and (1, 1)), one representing each class as per the GT.

Figure 7.

NFINDR extracted endmembers from Samson data.

Though this algorithm is mostly referred into the literature, there were no well-defined criteria for identifying the number of endmembers to extract. This challenge is tried to recover in [24] by using the notion of VD for defining the number of endmember available within the scene. Another issue with N-FINDR was its random initialization nature, which not only affects the algorithm convergence rate but also affects the final outcomes. There are several modifications are done in N-FINDR and made available in various versions [6, 25].

2.5 Simplex growing algorithm (SGA)

The simplex growing algorithm (SGA) [11] is a simplex based sequential pure pixel identification mechanism. It identifies the simplex having the maximum volume each time the new vertex is added. It is also called as the modified version of N-FINDR algorithm. The step-wise SGA is elaborated below:

2.5.1 Initialization

The SGA uses the VD estimation count to define the number of endmembers p to generate; and uses the e ₁ identified by the very first pure pixel selection procedure as the preferred initial pure pixel and initialize n = 1.

At $n\ge 1$ and for every data sample r, we compute the volume $V\left(e_1,\cdots ,e_n,r\right)$ defined by equation

$$V\left(e_1,\cdots ,e_n,r\right)=\frac{\left|det\left[\begin{array}{c}1\\ e_1\end{array}\begin{array}{c}1\\ e_2\end{array}\begin{array}{c}\cdots \\ \cdots \end{array}\begin{array}{c}1\\ e_n\end{array}\begin{array}{c}1\\ r\end{array}\right]\right|}{n!}$$ E10

It is the volume of simplex spanned by the vertices $e_1,e_2,\cdots ,e_n,r,$ symbolized with $S\left(e_1,e_2,\cdots ,e_n,r\right)$ . As the matrix $\left[\begin{array}{c}1\\ e_1\end{array}\begin{array}{c}1\\ e_2\end{array}\begin{array}{c}\cdots \\ \cdots \end{array}\begin{array}{c}1\\ e_n\end{array}\begin{array}{c}1\\ r\end{array}\right]$ in Eq. (10) is not essentially a square shaped matrix, the efficient DR technique like principal component analysis (PCA) or maximum noise fraction (MNF) is used to minimize the data dimensionality L to dimension n.

Identify the $e_{n+1}$ that produces the maximum of Eq. (10) and given by Eq. (11) given below:

$$e_{n+1}=arg\left\{\underset{r}{max}\left[V\left(e_{e1},\cdots ,e_n,r\right)\right]\right\}$$ E11

Stopping rule: if $n<p$ , then $n\leftarrow n+1$ and go to the step b. If not, the final set of vertices given by $\left\{e_1,e_2,\cdots ,e_p\right\}$ is considered as the p number of desired endmembers.

This is the step-wise process of SGA to find endmembers from hyperspectral images. The SGA algorithm efficiently solves the three significant challenges observed in N-FINDR algorithm. First, SGA performs the automated estimation of number of pure pixels to find, second, it produces the consistent final set of endmembers and third, it performs efficiently with less computational complexity. With the coincidence, the basic idea of using growing simplexes to find endmembers in SGA is analogous to the that used in VCA, but their corresponding approached are dissimilar. Most specifically, the VCA algorithm is unable to resolve the second challenge and produces inconsistent outcomes due to its random initialization nature [11].

2.6 Convex cone analysis (CCA)

The convex cone analysis (CCA) [12] algorithm assumes that physical quantities like radiance of reflectance of the hyperspectral images are always positive. The sample vectors produced by the discrete radiance or reflectance spectra can be articulated as linear mixtures or combinations of non-negative components, which reside within the nonnegative also called as convex region. The fundamental objective of CCA algorithm is to identify the boundary points for that convex region. The theoretical implementation of CCA algorithm is discussed below:

1. To practically implement this idea, the algorithm computes the eigenvectors of the data sample correlation matrix of the scene, and only chooses those eigenvectors conforming to the E largest eigenvalues (wherein, E is a previously defined number of endmembers to extract).

2. The algorithm then looks forward to identify the boundary points of the convex cone, wherein the linear combinations of these eigenvectors forms sample vectors that are exactingly non-negative, by using the equation below:

   $$h\left(x,y\right)=p_1+a_1{p}_2+\cdots +a_{E-1}{p}_E\ge 0$$ E12

   where $h\left(x,y\right)$ represents the spectral signature at pixel having spatial coordinates (x, y), the $p_i$ indicates the eigenvectors belonging to the largest eigenvalues, and 0 indicates the zero vector.

3. These identified points characterize the corner points of the convex cone region. Then, Eq. (6) is rewritten as:

   $$h\left(x,y\right)=\left[p_1\cdots p_E\right]\left[\begin{array}{c}1\\ a_1\\ \cdots \\ a_{E-1}\end{array}\right]=P_a\ge 0$$ E13

   wherein, the $p_i$ are the N-dimensional column vectors. For $N>E,P_a=0$ is an overdetermined system of the linear equations. If the elements of P are considered as the coefficients, and elements of a are considered variables, then, there are N equations of the given form:

   $$p_{j1}+a_1{p}_{j2}+\cdots +a_{E-1}{p}_{\mathit{jE}}=0_{,\text{for}}j=1,\cdots ,N$$ E14

   That defines (E-1)-dimensional hyperplane within the E-dimensional space.

4. Precise endmember set can be found by using the (E-1)-tuples from the N equations. These solutions will form the linear combinations of eigenvectors that is having minimum E-1 zeros.

5. At the end, the boundary points of the convex cone region are considered as the endmember that satisfy Eq. (12) or equivalently, Min[h (x,y)] = 0, where the smallest is taken over all $h_i\in h\left(x,y\right),i=1,\cdots ,N$ .

These sample vectors can be coined as final set of endmembers.

2.7 Endmember extraction using sparse component analysis (EESCA)

The endmember extraction using sparse component analysis (EESCA) [14] uses the sparse characteristics of abundance for endmember extraction from hyperspectral imagery. This algorithm is basically dependent of the primary endmember set of the VCA algorithm and thereafter two noteworthy iterative approaches are followed. The first one used to modify the endmember matrix and second one is to improve the pure pixel signature to enhance the algorithm accuracy. The step-wise mechanism of EESCA algorithm is summarized below:

The algorithm takes hyperspectral image as an input in the form of $X=\left[x\left(1\right),x\left(2\right),\dots ,x\left(T\right)\right]\in {\mathbb{R}}^{L\times T}$ and proceeds further.

1. Initiate the endmember matrix, i.e., $A^{\left[0\right]}$ using VCA algorithm.

2. Modify the pure pixel signatures using hyper-line estimation as per the steps given below:

   1. Set the value of error tolerance, i.e., $\epsilon$ , and repeat $k=1:K$

   2. Compute the distance $d\left(x\left(t\right),a_j^{\left[k-1\right]}\right)$ of all the data samples using Eq. (4) of [14], and allocate them in dissimilar classes $\left\{{\Omega }_j^{\left[k\right]},j=1,\cdots ,p\right\}$ .

   3. Choose a subset matrix of M pixels from ${\Omega }_j^{\left[k\right]}$ . Now, use SVD to modify the value of $a_i^{\left[k\right]}$ and eigenvalue ${\lambda }_i^{\left[k\right]}$ .

   4. If $\left|{\lambda }_j^k-{\lambda }_j^{\left[k-1\right]}\right|\le \epsilon$ , the pure pixel signature $a_j^{\left[k\right]}$ will be stable and not modified in the upcoming iteration.

3. Enhance the pure pixel signature using K-SVD as per the steps given below:

   1. State the matrix ${\omega }_p=\left\{i\left|1\le i\le N,s_p\left(i\right)\ne 0\right.\right\}$ to note down the non-zero indexes in s_p and define the maximum iteration time $t_{\max }$ .

   2. Compute E_p , the complete representation error matrix.

   3. Confine E_p by selecting the columns vectors conforming to ${\omega }_p$ , and derive the $E_p^R$ , the restricted error matrix.

   4. Pick he modified pure pixel column $a_p$ and update the corresponding coefficients vector using equation by Eq. (15) till the value of iterations reaches to t _max.

      $${‖X-\mathit{AS}‖}_F^2={‖X-\sum _{j=1}^p{a}_j{s}_T^j‖}_F^2={‖\left(X-\sum _{j\ne p}^p{a}_j{s}_T^j\right)-a_p{s}_T^p‖}_F^2={‖E_p-a_p{s}_T^p‖}_F^2$$ E15

At the end, the endmember matrix, i.e., $A=\left[a_1,a_2,\cdots ,a_p\right]\in R^{L\times P}$ , is given as final outcome. As per the comparative evaluation done in [14], it is found that EESCA is more capable and robust to deal with the mixed data and noise respectively. Along with its advantages, it has also few downsides. EESCA neither considers the spatial features nor touch the number of endmember estimation problem.

In further portion, three SI based EE approaches (i.e., ACO, PSO and ABC) are discussed. Though these algorithms have capability to solve the problem of combinatorial optimization, their common drawback is that they are computationally expensive due to the random searching mechanism [26, 27].

2.8 Ant colony optimization (ACO)

The ant colony optimization (ACO) [15] was first employed for endmember extraction from hyperspectral imagery in 2011 and termed as ACOEE. In ACOEE the problem of decomposition of mixed pixels is transformed into a problem of optimization and forms a feasible solution space to estimate the practical implication of the objective function. To generate the solution space and the heuristic information for endmember extraction a directed and weighted graph G is constructed, where G corresponds to each pixel in the HSI data. The route generated by artificial ant (i.e., feasible solution) contains m unlike vertices in G, where, m represents the number of endmembers. On the arrival of artificial ant at the vertex $v_i$ after $\left(t-1\right)$ times of moving, the probability of moving from $v_i$ to vertex $v_j$ is defined as [15],

$$p_{\mathit{ij}}^k\left(t\right)=\frac{{\tau }_{\mathit{ij}}^{\alpha }{\eta }_{\mathit{ij}}^{\beta }}{\sum _{j\in {\mathit{\text{allowed}}}_t}{\tau }_{\mathit{ij}}^{\alpha }{\eta }_{\mathit{ij}}^{\beta }},\forall j\in {\mathit{\text{allowed}}}_t$$ E16

where ${\tau }_{\mathit{ij}}$ shows the amount of pheromones in the edge $⟨v_i,v_j⟩$ and ${\mathit{\text{allowed}}}_t$ indicates the set of pixels covered by ant from $v_i$ (i.e., all pixels in the image except pixels that have been travelled by the ant) at time t. Parameters $\alpha$ and $\beta$ shows relative significance of pheromones and visibilities in the selection of route respectively, where the pheromones concentration is initialized with similar value. On the basis of path travelled by an artificial ant (i.e., endmembers), a remixed image can be built by using abundance estimation. The root-mean-square error (RMSE) between the original and the remixed image is used as the objective function to calculate the endmembers set. In the kth iterative sequence, $n$ ants creates dissimilar paths and pheromones in the edges updated using following equation,

$${\tau }_{\mathit{ij}}^{k+1}={\mathit{\rho \tau }}_{\mathit{ij}}^k+\Delta {\tau }_{\mathit{ij}}^k$$ E17

where $\Delta {\tau }_{\mathit{ij}}^k$ is the pheromone increment and $\rho$ is the factor of pheromone dissipation. If the smallest RMSE value in given repetition is $f_k$ and its equivalent path is ${\mathit{\text{path}}}_k$ , then $\Delta {\tau }_{\mathit{ij}}^k$ is,

$$\Delta {\tau }_{\mathit{ij}}^k=\left\{{}_{0⟨v_i,v_j⟩\notin {\mathit{\text{path}}}_k}^{Q/f_k⟨v_i,v_j⟩\in {\mathit{\text{path}}}_k}\right.,$$ E18

where Q represents constant which controls the rate of change of pheromones $\Delta {\tau }_{\mathit{ij}}^k$ , within acceptable range. The algorithm breaks when it reaches at the similar optimal path in multiple serial repetitive sequences or to the maximum number [15, 27].

ACO EE outperforms the popular VCA and NFINDR EE algorithm and improves the endmember representation [15].

2.9 Particle swarm optimization (PSO)

The particle swarm optimization (PSO) [16] technique is a mostly used swarm intelligence technique to handle the problem of global optimization. Discrete PSO for endmember extraction (DPSO EE) was firstly adapted by Zhang et al. in 2011, where concept of PSO was utilized. In DPSO EE the combinatorial optimization problem is resolved by converting endmember extraction into feasible solution space and objective function, i.e., PSO searches for endmembers in the distinct feasible solution by outlining the particles locations and their velocities along with binary coding, and shows better endmember extraction results for HSI data than state of art methods [16, 28].

Linear spectral mixture model was used for spectral unmixing as shown in Eq. (19),

$$r_i=\sum _{j=1}^m{e}_j{\alpha }_{\mathit{ij}}+{\epsilon }_i$$ E19

where ${\left\{r_i\right\}}_{i=1}^n$ represents L bands and n pixels remote sensing image with $r_i$ as column vector of spectrum of ith pixel, ${\left\{e_i\right\}}_{i=1}^n$ shows endmember set, ${\epsilon }_i$ is the random error and ${\alpha }_{\mathit{ij}}$ is the abundance of jth endmember in ith pixels. The least square method has been utilized for calculating area ratio of each endmember in a pixels, i.e., ${\alpha }_{\mathit{ij}}$ .

If the values of ${\left\{r_i\right\}}_{i=1}^n$ and ${\left\{e_i\right\}}_{i=1}^n$ are known then the image can remixed as, $r_i=\sum _{j=1}^m{e}_j{\alpha }_{\mathit{ij}},$ $i=1,2,\dots ,n$ . The root mean square error (rmse) between an original image and remixed image is calculated as [29],

$$\mathit{\text{rmse}}\left({\left\{r_i\right\}}_{i=1}^n,{\left\{e_i\right\}}_{i=1}^n\right)=\frac{1}{n}\sum _{i=1}^n\sqrt{\frac{1}{L}}{‖r_i-{\stackrel{⌢}{r}}_i‖}_2^2$$ E20

If the value of rmse is less, then the result of endmember extraction is better. Therefore, endmember extraction can be defined as the combinatorial optimization problem as,

$$\begin{array}{l}Min\mathit{\text{rmse}}\left({\left\{r_i\right\}}_{i=1}^n,E\right)\\ \mathrm{s}.\mathrm{t}.E\in C\left({\left\{r_i\right\}}_{i=1}^n,E,m\right)\end{array}$$ E21

Here, $C\left({\left\{r_i\right\}}_{i=1}^n,E,m\right)$ referred as, set of subsets of ${\left\{r_i\right\}}_{i=1}^n$ that contains m elements. $\mathit{\text{rmse}}\left({\left\{r_i\right\}}_{i=1}^n,E\right)$ is the objective function and $C\left({\left\{r_i\right\}}_{i=1}^n,E,m\right)$ is the feasible solution space.

The feasible solution space of the optimization problem is used to search particles by PSO. The adaptability function is obtained on the basis of constantly moving particles in the feasible space. In case of discrete feasible solution space D-PSO is improved to search particles in it on the basis of PSO. Moreover, the mapping relationship in the image and the feasible solution space can be given as [29],

$$\begin{array}{l}G:C\left({\left\{r_i\right\}}_{i=1}^n,m\right)\to X_{n,m}\\ E\mapsto \left(x_1,x_2,\dots x_n\right).\end{array}$$ E22

where $X_{n,m}=\left\{\left(x_1,x_2,\dots x_n\right)|x_i\in \left\{0,1\right\},\sum _{i=1}^n{x}_i=m\right\}$

From above formula, if a pixel $r_i$ in ${\left\{r_i\right\}}_{i=1}^n$ is selected as an endmember, then value of x_i corresponding to r_i in x is 1; otherwise, the value is 0.

The DPSO EE shows better performance compared to VCA and NFINDR, but it has few limitations. The limitation includes several additional parameters that affect the performance of optimization and no assurance during about whether the location of selected particle is best or not [16].

2.10 Artificial bee colony (ABC)

The artificial bee colony (ABC) [17] algorithm is utilized for solving optimization problems by using the searching activity of bees in nature. The colony’s search space is considered as feasible solution space for problem which is to be solved, where feasible solution is referred as a food source and quantity of nectar present in each source is referred as its fitness, which is linked with objective function value created by corresponding solution. Bees can be categorized into three types: employed bees, scout bees and onlooker bees. Xu Sun et al. used artificial bee colony for endmember extraction in 2015 [17, 27]. The step-wise workflow of ABC algorithm is given below:

1. The position of every employed bee is considered as position of a food source. The data about particular food source is recorded and fitness of respective food source is calculated. After finding with adjacent food sources, if fitness of newly found food source is better than previous one then update the former food source with new one otherwise continue to search.

2. If $X_i={\left(X_{i1},X_{i2},\dots ,X_{\mathit{iM}}\right)}^T$ represents the ith source food (i.e., location of employed bee) then adjacent searching can be defined as follows:

   $${\phi }_{\mathit{ij}}=x_{\mathit{ij}}+\Delta \left(x_{\mathit{ij}}-x_{\mathit{sj}}\right)$$ E23

   where, ∆ shows a random value in the range of [−1, 1], $s\ne i,s\in \left\{1,2,\dots M\right\}$ and $j\in \left\{1,2,\mathrm{..},M\right\}$ the greedy selection operator selects a food source with a better fitness value.

3. Each onlooker bee selects a food source on the basis of fitness acquired by the employed bees. It repeats similar policy of employed bee. The “follow probability” of an onlooker bee of selecting the jth food source is calculated as follows:

   $$p_j=\frac{\mathit{fit}\left(x_j\right)}{{\sum }_{i=1}^{\mathit{num}}\mathit{fit}\left(x_i\right)}$$ E24

   where $\mathit{num}$ shows the number of food sources which is greater than M and $\mathit{fit}\left(x_i\right)$ shows the fitness function of the ith food source and.

The fitness of a food source is determined using Eq. (25), where $f\left(x_i\right)$ is the objective function of the given optimization problem:

$$\mathit{fit}\left(x_i\right)=\frac{1}{f\left(x_i\right)}$$ E25

The scout bee randomly reaches at a food source in the respective feasible solution space, converting itself into an employed bee to calculate the fitness value of the respective food source. Each random search can compute the best fitness position using the iterative process. If a food source has not been updated in a long time, then it is considered a global optimal source [17]. The ABC algorithm needs to face higher computational complexity for the images of larger sizes containing more number of endmembers.

3. Conclusion

The ultimate objective of using HSI data is to achieve the higher material classification or mapping accuracy, but it is adversely affected due the ground truth unavailability or spectral mixing problems. To overcome the same efficiently and enhance the HSI data classification or mapping accuracy, an image derived endmembers offers a precise and useful solution. To extract the endmembers from HSI data, several approaches based on either convex geometry or statistical information were tried over the past few years. Out of which, due to the higher volume of the data, the convex geometry based EE algorithms were found more effective and widely accepted and discussed in the literature. The popular algorithms comes under this criteria are PPI, VCA, SMACC, NFINDR, CCA and SGA. Apart from these techniques, some other statistical information and artificial intelligence based approaches are also well available in the literature, but unfortunately, none of the method provides the automated, data and application independent outcomes. Every method has its own limitations and challenges. These methods provides either inconsistent or misleading outcomes due to the random algorithm initialization nature, unclear criterion to define the values of additional input parameters and weakly defined algorithm stopping rules. In case of PPI, the algorithm is very sensitive to the value of iterations and threshold, where higher number of iterations results into higher computational complexity and unclear criteria to set the threshold value leads to the endmember selection problem. In case of VCA, SMACC and NFINDR, the basic requirement is that, the number of unique materials available into the image needs to be known primarily, which is hardly possible. To satisfy this requirement, some EE algorithms have opted the outcome of VD estimation (which gives the number of distinct signal sources available within the scene) as a reference to extract the number of endmembers. Further investigation of VD estimation techniques revealed that, these methods are also produce misleading results. More detailed and comparative analysis between various endmember extraction algorithms is performed time to time in [30, 31, 32, 33]. Therefore, to overcome most of the challenges available in recent EE algorithm, there is a great need of a fully automated and robust endmember extraction algorithm. The efficient endmember extraction algorithm need to be data and application independent and must be able to find precise set of endmembers without any prior knowledge.

Acknowledgments

This work has been supported by DST, GOI, [under MRP No. BDID/01/23/2014-HSRS/35 (ALG-V)], UGS SAP II DRS Phase II, DST FIST and Department of Computer Science and Information technology, Dr. Babasaheb Ambedkar Marathwada University, Aurangabad-(MS), India.