Student Member, IEEE, Julien Mairal, Senior Member, IEEE Alexandre Zouaoui

Fellow, IEEE Jocelyn Chanussot

Julien Mairal

Image Processing and Machine Learning for Hyperspectral Unmixing: An Overview and the HySUPP Python Package

Keywords: Hyperspectral, unmixing, endmember extraction, abundance estimation, linear mixture, machine learning, deep learning, optimization

Spectral pixels are often a mixture of the pure spectra of the materials, called endmembers, due to the low spatial resolution of hyperspectral sensors, double scattering, and intimate mixtures of materials in the scenes. Unmixing estimates the fractional abundances of the endmembers within the pixel. Depending on the prior knowledge of endmembers, linear unmixing can be divided into three main groups: supervised, semisupervised, and unsupervised (blind) linear unmixing. Advances in image processing and machine learning substantially affected unmixing. This paper provides an overview of advanced and conventional unmixing approaches. Additionally, we draw a critical comparison between advanced and conventional techniques from the three categories. We compare the performance of the unmixing techniques on three simulated and one real dataset. The experimental results reveal the advantages of different unmixing categories for different unmixing scenarios. Moreover, we provide an open-source Python-based package available at https://github.com/BehnoodRasti/HySUPP to reproduce the results.

I. INTRODUCTION

S PECTRAL unmixing is a crucial processing technique in hyperspectral remote sensing that can play a vital role in various fields such as mineral exploration, agriculture and crop monitoring, environmental monitoring, urban planning, remote sensing of planetary surfaces, pollution monitoring, medical imaging, water quality assessment, etc. The ability to separate and identify different materials in an image is made possible by the contiguous spectra captured by hyperspectral sensors. Through the use of endmembers, which are the unique spectral signatures of materials, unmixing algorithms can decompose the mixed spectral data into its constituent parts. However, due to low spatial resolution, multiple scattering, and intimate mixing, the measured spectrum within a pixel is generally a complex mixture of the pure spectra of the constituent materials, making unmixing a challenging task. Fig. 1 demonstrates how the reflectance of a mixed pixel captured by an optical hyperspectral camera is composed of two endmembers within that pixel. In hyperspectral remote sensing, a mixing model represents the observed spectral pixel as a function of the endmembers and their corresponding fractional abundances within the pixel's area. Unmixing is the process of estimating the fractional abundances, either by estimating or extracting the endmembers or by relying on a library of endmembers. It may also involve determining the number of endmembers present. The mixing model is either linear or nonlinear, depending on the interaction of the incident light and the materials in the scene or sample.

In linear unmixing, the endmembers are assumed to be linearly mixed, which is valid when each light ray interacts with only one material before reaching the sensor, as shown in Figure 2 (a). This assumption is common in Earth observation applications where macroscopic problems exist. The spatial resolution of the sensor plays a crucial role in macroscopic scenarios, where a pixel may contain multiple materials. In such cases, the pixel is a mixture of different materials, such as tree and soil, as shown in Figure 2 (a).

Another assumption is bilinear mixing, which assumes double scattering or that the light ray interacts with two materials before reaching the sensor, as shown in Figure 2 (b). However, in microscopic scenarios where the pure materials are intimately mixed within a pixel, and the light undergoes multiple scattering and reflections by several materials, the linear approximation fails, and nonlinear models that are more complex than bilinear should be used [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], [START_REF] Dobigeon | Nonlinear unmixing of hyperspectral images: Models and algorithms[END_REF].

In macroscopic scenarios, nonlinearity is also induced by noise, atmospheric effects, temporal effects, illumination variations, and materials' intrinsic variation, which may cause spectral variability shown in Fig. 3 [START_REF] Borsoi | Spectral variability in hyperspectral data unmixing: A comprehensive review[END_REF]. Illumination variations are mainly due to two effects: varying terrain topography, which affects the angles of the incident radiation, and the occlusion of the light source by other objects (leading to shaded areas). The intrinsic variation of spectral signatures of the materials is mainly due to physicochemical changes. For instance, soil's signature might change dramatically by variations in its composition and moisture content. Another example is the leaf's signature changes throughout the temporal seasons. Scale differences can cause spectral variability as well. For instance, a tree contains leaves, branches, fruits, bark, bark, fruits or flowers, etc., leading to the intrinsic variability of the tree spectra at different scales. The atmospheric effect should be corrected, and noise can be removed [START_REF] Coquelin | Hyde: The first open-source, python-based, gpu-accelerated hyperspectral denoising package[END_REF]- [START_REF] Rasti | Image Restoration for Remote Sensing: Overview and Toolbox[END_REF].

The unmixing process is subjective and depends on how the endmembers are defined. For example, when detecting buildings, the reflectance of a brick could be considered an endmember. However, bricks are composed of several materials, such as clay, sand, and concrete, making it challenging to define the endmembers of interest for detecting specific materials or ratios of different materials inside bricks.

In the former, linear unmixing may be sufficient to solve the problem. However, for intimate mixtures and the latter, nonlinear models must be used. For example, the Hapke model [START_REF] Hapke | Theory of Reflectance and Emittance Spectroscopy[END_REF], [START_REF] Hapke | Bidirectional reflectance spectroscopy: 1. theory[END_REF] suggests that mixing occurs at the albedo level rather than the reflectance level. Therefore, more complex nonlinear models are necessary for such cases to accurately unmix the spectral data. It is worth mentioning that, unlike the common assumption, the proportion of an endmember within a pixel or in a scene is not the percentage of the endmember material within that pixel or in that scene [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF]. According to Hapke [START_REF] Hapke | Theory of Reflectance and Emittance Spectroscopy[END_REF], [START_REF] Hapke | Bidirectional reflectance spectroscopy: 1. theory[END_REF], the proportion of an endmember in a linear mixture indicates the relative area of that endmember. Therefore, reflectance is not a linear mixture of the mass or cross-sectional area of the endmember materials. Despite this, linear unmixing has shown significant value in remote sensing applications in the past decades [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF]. However, when associating an unmixing problem with an application, we should be aware of the simplifications of the model we use.

Assuming three materials in a scene correspond to three endmembers (e 1 , e 2 , and e 3 ), the observed data can be projected into a 2D subspace where all the data points are enclosed in the data simplex formed by the convex hull of endmembers (see Fig. 11). In contrast, in a nonlinear scenario, the simplex does not contain all the data points. This paper mainly considers linear EO applications and macroscopic problems and focuses on linear mixing models. However, we briefly discuss the bilinear and nonlinear unmixing approaches.

A. Key trends in the hyperspectral unmixing realm

Hyperspectral unmixing has undergone various developmental stages, as illustrated in figure 4. In its early years, from 1992 to 2008, this research topic received minimal attention, with an annual publication count of no more than 18. However, starting from 2009 and continuing until 2016, there was a steady increase in interest, with the number of publications even surpassing 200 in 2016. Since 2017, hyperspectral unmixing has remained a prominent and dynamic area of research, driven by the emergence of deep learning techniques for unmixing. In 2022, this field reached its fourth most prolific year on record, with 164 publications.

During the period when hyperspectral unmixing gained popularity, the primary scientific programming language of choice was MATLAB, as depicted in figure 5. However, in recent times, there has been a steady increase in interest and adoption of Python, while MATLAB's popularity has waned. This shift can be attributed, in part, to the growing trend of open science, where open-sourcing scientific code has become widely embraced. Consequently, the present scenario suggests that developing and releasing an open-source unmixing package in Python is more advantageous than doing so in MATLAB, even though it may involve the translation of existing unmixing code from MATLAB to Python.

B. Existing Overview Publications and Unmixing Packages

An early survey on hyperspectral unmixing was given in [START_REF] Keshava | A Survey of Spectral Unmixing Algorithms[END_REF], which discusses basic geometrical and statistical methods. In [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], linear model-based unmixing techniques were divided into three categories: geometrical, statistical, and sparse regression-based approaches. A Matlab toolbox is available at https://openremotesensing.net/. However, the toolbox is incomplete, and some methods, such as dependent component analysis (DECA), are missing. In recent years, deep learning and neural networks have become state-of-the-art for many tasks in machine learning and image processing. Consequently, many unmixing approaches were proposed based on shallow and deep neural networks. A comparison of autoencoderbased networks was drawn in [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF]. The authors discussed autoencoder-based architectures divided into five categories, i.e., Sparse Nonnegative Autoencoders, Variational Autoencoders, Adversarial Autoencoders, Denoising Autoencoders, and Convolutional Autoencoders. They further discuss the choice of different modules, such as different activation functions or loss functions, and they compared shallow networks to deep ones. They also provide a TensorFlow-based Python package that is available on GitHub. However, the package is limited to blind unmixing approaches based on autoencoders. It does not discuss or compare the supervised, semisupervised, and more conventional blind unmixing methods.

In [START_REF] Chen | Integration of physics-based and data-driven models for hyperspectral image unmixing: A summary of current methods[END_REF], some model-based and neural network-based unmixing approaches were explained but without experimental comparisons. A list of resources for the approaches was given.

A survey on endmember variability in Spectral Mixture Analysis (SMA) was given in [START_REF] Ben | Endmember variability in spectral mixture analysis: A review[END_REF]. In [START_REF] Zare | Endmember variability in hyperspectral analysis: Addressing spectral variability during spectral unmixing[END_REF], an overview of unmixing methods that address endmember variability was given. A comprehensive overview of the unmixing methods that address spectral variability was recently provided in [START_REF] Borsoi | Spectral variability in hyperspectral data unmixing: A comprehensive review[END_REF], and a list of Matlab codes was also given. In [START_REF] Plaza | Recent Developments in Endmember Extraction and Spectral Unmixing[END_REF], an overview of endmember extraction approaches was given. Review papers on hyperspectral remote sensing data analysis briefly discussed the unmixing methods [START_REF] Bioucas-Dias | Hyperspectral remote sensing data analysis and future challenges[END_REF], [START_REF] Ghamisi | Advances in hyperspectral image and signal processing: A comprehensive overview of the state of the art[END_REF]. In [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF], a survey of nonlinear unmixing methods is given.

There are other open-source tools such as HyperMix [START_REF] Jimenez | Hypermix: An open source tool for hyperspectral imaging[END_REF], Spectral Python (SPy) (https://www.spectralpython.net/), Spectral Library Tool (https://spectral-libraries.readthedocs.io/en/latest/), PySptools (https://pypi.org/project/pysptools/), that include basic algorithms for estimation of the number of endmembers, endmember extraction, abundance estimation, and some library tools, and library-based methods. Thus, there is a need for a comprehensive package that covers the methodologies across different unmixing categories and contains state-ofthe-art image processing and machine learning techniques.

C. Significance of Contribution

This paper aims at the followings:

• Providing researchers with a comprehensive yet technical overview of all the essential topics regarding linear unmixing techniques. • Categorizing the unmixing approaches considering prior knowledge available about endmembers. Linear unmixing can be divided into three main categories: supervised, semi-supervised (library-based), and unsupervised (blind) unmixing.

• Compare the unmixing methods in terms of prior knowledge of the endmembers and draw conclusions that can help researchers to select an appropriate unmixing method to tackle real-world challenges. We compare conventional and deep learning-based unmixing approaches in those categories for three simulated and two realworld datasets. For the simulated datasets, we consider three scenarios: a simple, pure pixel dataset, a dataset with spectral variations, and a challenging dataset with no pure pixel. Such comparisons provide insight to researchers into which category to use for their application. Additionally, the comparisons reveal the drawbacks of the categories, which motivate the developers to investigate new ideas to address them.

•

We provide an open-source HyperSpectral Unmixing Python Package (HySUPP). HySUPP is the first opensource python-based hyperspectral unmixing package to include supervised, semi-supervised, and blind unmixing methods. The package will benefit the geoscience and remote sensing community, including researchers, developers, lecturers, and students. The package installation is straightforward since HySUPP relies on a few dependencies. In addition, all the methods can be run using a single command line instruction.

II. HYSUPP: HYPERSPECTRAL UNMIXING PYTHON PACKAGE

HySUPP exhibits a list of highly desirable properties summarized as follows: i) completeness, ii) reproducibility, iii) extensibility and iv) homogeneity. It implements common best practices and enables simple benchmarking of unmixing techniques thanks to user-friendly command line instructions.

A. Features a) Completeness: As a practitioner, the ability to explore and experiment with different unmixing techniques is crucial since no single approach can consistently outperform others in all unmixing scenarios. Thus, ensuring the completeness of our toolbox becomes essential. Our toolbox is designed to cover all three types of unmixing -supervised, semi-supervised, and unsupervised -while striving to be as representative as possible of the various unmixing approaches, although aiming for exhaustive inclusion would be impractical. Currently, HySUPP provides access to a diverse set of 20 different unmixing methods, including 6 supervised, 6 semi-supervised, and 8 unsupervised techniques.

b) Reproducibility: Ensuring experimental reproducibility is crucial when exploring various unmixing techniques, as it guarantees the robustness of the conclusions drawn by users. In line with this principle, our toolbox offers the ability to seed experiments, facilitating reproducibility through repeatable noise generation. Additionally, HySUPP automatically saves estimates outputs, providing users with a convenient way to review and compare results, thus enhancing the reliability of their research findings. c) Extensibility: HySUPP's architecture is designed to support the effortless integration of new methods, providing a platform for future advancements in hyperspectral unmixing. Leveraging configuration files, we empower researchers to conduct experiments with ease and flexibility, enabling them to effortlessly incorporate their own models into the toolbox.

d) Homogeneity: The uniformity of HySUPP's codebase ensures consistency across different components, enhancing its overall usability and reliability. More specifically, we use common methods for models and common attributes for datasets.

e) Best practices: We offer a straightforward yet potent method to monitor the objective function of each approach using Python's tqdm library. Furthermore, our pipeline incorporates precise endmember auto-alignment through the utilization of the munkres algorithm, provided ground truth abundance maps are accessible. This enhancement ensures the accurate computation of unmixing performance. By establishing these best practices, we streamline the implementation process and foster a collaborative environment where researchers can easily build upon existing work and share their contributions effectively.

f) Benchmarking: Owing to its rich set of features, HySUPP enables users to easily benchmark all implemented methods on their dataset of choice. Our toolbox currently provides 3 synthetic datasets corresponding to different unmixing scenarios. Moreover, we incorporate 4 distinct metrics to thoroughly evaluate unmixing accuracy across methods. Finally, leveraging mlxp [START_REF] Arbel | Mlxp: A framework for conducting machine learning experiments in python[END_REF] results query tool, we empower users to analyze and visualize their results in an appealing and informative manner.

B. Example

The following command line instruction provides an example on how to run a semi-supervised unmixing technique, SUnCNN [START_REF] Rasti | SUnCNN: Sparse unmixing using unsupervised convolutional neural network[END_REF], on the DC1 dataset using an optional custom value for the signal-to-noise ratio (SNR):

$ python unmixing.py mode=semi data=DC1 model=SUnCNN noise.SNR=30

Table I lists unmixing methods included in HySUPP with their corresponding dependencies. We highlighted our main contributions and the link to the original implementations.

III. LINEAR UNMIXING

For frequent Earth observation, hyperspectral sensors are installed on satellites. The light passes through the atmosphere before reaching the target and then passes through the atmosphere to reach the sensor. The atmosphere absorbs and scatters light, and therefore, atmospheric corrections should be applied in addition to converting radiance to reflectance. Arguably, atmospheric corrections are unnecessary for airborne, drone-borne, and close-range imaging.

Fig. 6 (a) simplifies the concept of sensing using a satellite hyperspectral sensor. As can be seen, the sensor captures a pixel that contains three materials, i.e., Water, Tree, and Soil. At the sensor, irradiance is corrected to the reflectance in Fig. 6 (b). Note that this reflectance is between zero and one. Fig. 6 (c) shows the concept of linear unmixing, which demonstrates that the measured pixel contains 20% water (blue endmember), 50% soil (yellow endmember), and 30% tree (green endmember). Arguably, these portions cannot be negative, known as abundance nonnegativity constraint (ANC), and the summation of those portions is 100%, known as abundance sum to one constraint (ASC). Additionally, the endmembers are reflectance and cannot be negative (indeed, the endmembers have to be bounded between 0 and 1), which led to the endmember nonnegativity constraint. Here, we discuss two widely used linear mixing models, i.e., the lowrank linear mixing model (so-called LMM) and the sparse and redundant linear model. The former was mainly used for supervised and blind unmixing and the latter for sparse unmixing. Some other linear models have been used in the Refactor into a SupervisedUnmixingModel sub-class SiVM [START_REF] Heylen | Fully constrained least squares spectral unmixing by simplex projection[END_REF] github numpy Refactor into a BaseExtractor sub-class SISAL [START_REF] Bioucas-Dias | A variable splitting augmented lagrangian approach to linear spectral unmixing[END_REF] github numpy Refactor into a BaseExtractor sub-class UnDIP [START_REF] Rasti | UnDIP: Hyperspectral unmixing using deep image prior[END_REF] github torch Refactor separate scripts into a single model VCA [START_REF] Nascimento | Vertex component analysis: A˜fast algorithm to extract endmembers spectra from hyperspectral data[END_REF] github numpy Replace scipy dependency by numpy. Add random seed CLSUnSAL [START_REF] Iordache | Collaborative sparse regression for hyperspectral unmixing[END_REF] github numpy Implementation based on SUnSAL MUA SLIC [START_REF] Borsoi | A fast multiscale spatial regularization for sparse hyperspectral unmixing[END_REF] github numpy, skimage MATLAB code translation using skimage's SLIC S 2 WSU [START_REF] Zhang | Spectral-spatial weighted sparse regression for hyperspectral image unmixing[END_REF] github numpy, scipy MATLAB code translation SUnAA [START_REF] Rasti | Sunaa: Sparse unmixing using archetypal analysis[END_REF] github numpy, spams Refactor to match SemiUnmixingModel sub-class SUnCNN [START_REF] Rasti | SUnCNN: Sparse unmixing using unsupervised convolutional neural network[END_REF] github torch Refactor separate scripts into a single model SUnSAL [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF] github numpy Refactor into SemiUnmixingModel sub-class ADMMNet [START_REF] Zhou | Admm-based hyperspectral unmixing networks for abundance and endmember estimation[END_REF] torch Implemented from scratch BayesianSMA [START_REF] Dobigeon | Joint bayesian endmember extraction and linear unmixing for hyperspectral imagery[END_REF] webpage matlab.engine numpy Python wrapper around existing MATLAB code CNNAEU [START_REF] Palsson | Convolutional autoencoder for spectral-spatial hyperspectral unmixing[END_REF] github torch Convert existing keras implementation into torch EDAA [START_REF] Zouaoui | Entropic descent archetypal analysis for blind hyperspectral unmixing[END_REF] github torch Refactor to match BlindUnmixingModel sub-class MiSiCNet [START_REF] Rasti | Misicnet: Minimum simplex convolutional network for deep hyperspectral unmixing[END_REF] github torch Refactor separate scripts into a single model MSNet [START_REF] Yu | Multi-stage convolutional autoencoder network for hyperspectral unmixing[END_REF] github torch Refactor to match BlindUnmixingModel sub-class NMFQMV [START_REF] Zhuang | Regularization parameter selection in minimum volume hyperspectral unmixing[END_REF] github matlab.engine numpy Python wrapper around existing MATLAB code PGMSU [START_REF] Shi | Probabilistic generative model for hyperspectral unmixing accounting for endmember variability[END_REF] github torch Refactor to match BlindUnmixingModel sub-class literature, which will be discussed throughout the paper. We also discuss the components of a basic neural network model.

Unmixing techniques can be categorized into three main groups considering the prior knowledge of endmembers. Supervised and unsupervised (blind) unmixing use the same lowrank mixing model, while semisupervised unmixing uses the sparse and redundant linear model. The endmember library is overcomplete; consequently, the abundances are desired to be sparse. Blind unmixing estimates both endmembers and abundances simultaneously. Fig. 7 graphically compares those three groups. It is worth mentioning that unsupervised (blind) unmixing can be applied to at sensor irradiance, however, the data cannot be interpreted and the endmembers extracted/estimated cannot be associated with the corresponding materials. One can apply supervised unmixing to irradiance if the pure pixels are also irradiance. Semisupervised unmixing cannot be applied to sensor radiance since they rely on a spectral library.

A. Low-rank Linear Mixture Model

Every p-dimensional pixel y (i.e., the sensor has p bands) captured by a hyperspectral sensor can be represented as a linear combination of the endmembers within the pixel. Let us assume that matrix E ∈ R p×r contains r endmembers within the pixel. Then, y ∈ R p is represented as

y = Ea + n, s.t. r i=1 a i = 1, a i ≥ 0, i = 1, 2, .., r, (1)
where n denotes the p-dimensional random vector denoting the additive random Gaussian noise. Using matrix notations, we can represent all the pixels Y as

Y = EA + N, s.t. A ≥ 0, 1 T r A = 1 T n , (2)
where Y∈ R p×n is the observed HSI, with n pixels and p bands, N ∈ R p×n is noise, and A ∈ R r×n contain the r endmembers and their fractional abundances, respectively. 1 n indicates an n-component column vector of ones. Model 2 is known as the linear mixture model (LMM).

B. Low-rank LMM and Simplex Volume

Suppose that the endmembers are affinely independent i.e., e 2 -e 1 , ..., e r -e 1 are linearly independent, then

S {Ea ∈ R p | r i=1 a i = 1, a i ≥ 0, i = 1, 2, .., r} (3)
is (r -1)-Simplex. Indeed, this is the convex hull of the vertices, i.e., e i 's assuming no noise. Therefore, if we ignore the noise, any point of the dataset belongs to S (see Fig. 11 (a)). We should note that columns of A belong to (r -1)probability simplex

∆ r {a ∈ R r | r i=1 a i = 1, a i ≥ 0, i = 1, 2, .., r}. (4)
Therefore, a i ∈ ∆ r means ASC and ANC. Indeed, the vertices of ∆ r are r unit vectors in R r . Hereafter, we refer to (3) as data simplex and to (4) as abundance simplex. It is worth mentioning that, unlike data simplex, the abundance simplex is also valid for nonlinear models.

Let us assume S r is an r-simplex in R r i.e., e i ∈ R r for i = 0, 1, 2, .., r. Then the volume of this simplex is given by

V ol(S r ) = 1 r! |det[e 1 -e 0 , ..., e r -e 0 ]| = 1 r! |det 1 ... 1 e 0 ... e r |, (5)
where det is the determinant of a matrix. In the linear mixture model (LMM), we have an (r -1)-Simplex in R r (after the projection into an r dimensional subspace, see IV-A) i.e., e i ∈ R r for i = 1, 2, .., r (E r ∈ R r×r ). The volume of an (r -1)simplex is zero in R r . Therefore, the extended simplex, which contains the origin (i.e., E 0 = [0, E]), can be used to calculate the volume. Hence,

V ol(E 0 ) = 1 r! |det[e 1 , ..., e r ]| = 1 r! |det(E)|. ( 6
)
Alternatively, the data points can be shifted to the origin, and the volume of the (r -1)-simplex can be computed in R r-1 given by

V ol(E r ) = 1 (r -1)! |det[e 2 -e 1 , ..., e r -e 1 ]| = 1 (r -1)! |det 1 ... 1 e 1 ... e r |. (7)
It is worth mentioning that the projection into a subspace (IV-A) is necessary to form a squared matrix and calculate the determinant.

C. Sparse and Redundant Linear Mixture Model

The observed spectra can be represented using a sparse and redundant linear mixture model given by

Y = DX + N, s.t. X ≥ 0, 1 T m X = 1 T n , (8)
where D ∈ R p×m (p m) denotes the spectral library containing m endmembers and X ∈ R m×n is the unknown fractional abundances to estimate. Note that D is an overcomplete dictionary and therefore, should be a well-designed dictionary. A well-designed dictionary contains endmembers of the materials in the scene and can sparsify the redundant X. Hence, a spectral library can be pruned based on the spectral angles of the spectra (i.e., spectra with small degree differences will be removed). However, this comes with a risk of losing material endmembers if they are scaled versions of each other. In the case of a well-designed dictionary, the pixels are a mixture of a few atoms of the dictionary, and therefore X is a sparse matrix. Note that we get a zero row if an endmember material does not exist in the observed spectra. The rest of the rows will contain zeros since the abundances are often sparse. This model is often used in sparse unmixing. The fractional abundances X are estimated using sparsity-enforcing penalties/constraints in a sparse regression formulation.

D. Shallow/Deep Neural Network Model

Shallow or deep neural networks can be used for supervised, semi-supervised, or blind unmixing. The network is often based on autoencoder-based architecture, as shown in Fig. 8. The encoder E , which could be deep or shallow, encodes the spectral pixels into the abundances given by

a = E (y), (9)
and the decoder D reconstruct the pixel given by ŷ = D(a).

In every layer, an activation function is used to apply nonlinearity. The most common choices are Rectified Linear Unit ReLU (ReLU(x) = max(0, x)) and Leaky ReLU (LReLU(x) = max(ax, x), where a is very small, and a common choice is a = 0.1). Batch normalization (BN) is commonly used in each layer to speed up the learning process.

The BN function is given by

BN (x) = α 2 x -µ √ σ 2 + + α 1 , (11)
where α 2 and α 1 are learnable parameters, µ and σ are the mean and standard deviation of the batch, and is a minimal value.

A linear layer can be used for the decoder to reconstruct the pixel. The decoder weights are the endmembers fixed for supervised unmixing and can be learned through the learning process in the case of unsupervised (blind) unmixing. ASC and ANC can be realized differently; Applying nonnegative activation functions such as ReLU and normalization for ASC given by [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF] a = a

r i=1 a i , (12)
alternatively, one can use

a = |a| r i=1 |a i | , (13)
A common way is to use a softmax function sof tmax(a) = e a r i=1 e ai .

(

) 14
The scaled version of the softmax can be used as sof tmax(γa i ) where γ > 1 is constant. ASC and ANC can also be implemented by adding a penalty term to the loss function. However, it degrades the stability of the training process and decreases the convergence time, and therefore, it is not recommended [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF]. We should note that those functions are realized as the last hidden layer of the decoder. The loss function is often chosen as a reconstruction term plus a regularization term given by

L reg (y, ŷ) = L(y, ŷ) + λJ(a, Ŵen , Ŵde ), ( 15
)
where L is often selected as a reconstruction term such as 2 norm (or MSE)

L 2 (y, ŷ) = ||y -ŷ|| 2 2 , (16)
or

L SAD (y, ŷ) = arccos y T ŷ ||y|| 2 ||ŷ|| 2 , (17)
or

L SID (y, ŷ) = p i=1 p i log( p i q i ) + p i=1 q i log( q i p i ), (18)
where

p i = y i p i=1 y i , q i = ŷi p i=1 ŷi , (19)
or a combination of them. The regularizer J can be a combination of several functions, such as a Frobenius norm on W en with a sparsity norm such as 1 norm on a, and a Frobenius norm or a geometrical (MV induced) norm on W de . Note that SAD is scale-invariant and, therefore, can capture the spectral variability of the pixel. This is a big advantage for real-world datasets. On the other hand, utilizing SAD for datasets containing materials having scaled endmembers is not recommended since the network cannot be trained to distinguish those materials. Assuming a single-layer encoderdecoder, the abundances are given by

a = sof tmax(BN (RelU (W en y + b))), ( 20
)
where b denotes the biases of the layer. The reconstructed pixel is given by ŷ = W de a.

Fig. 8. A basic encoder-decode architecture used for unmixing. The (shallow or deep) encoder encodes the pixels to abundance. ASC and ANC can be enforced on the bottleneck. The shallow decoder reconstructs the pixels, and the decoder's weights are the endmembers.

E. The Effect of Noise in Unmixing

Hyperspectral data are usually downgraded by different sources, including atmospheric effects, imaging artifacts, and instrumental noise. These sources can further affect hyperspectral analysis, including unmixing [START_REF] Rasti | Noise reduction in hyperspectral imagery: Overview and application[END_REF]. Unmixing is an inverse problem that could be very sensitive to noise, mainly when endmembers are highly correlated. Therefore, unmixing techniques may fail due to noise. In [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], so-called signal-tonoise-ratio spectral distribution (SNR-SD) was suggested to determine if the unmixing inverse problem gives acceptable results. However, the hyperspectral denoising field has considerably evolved during the past decade [START_REF] Coquelin | Hyde: The first open-source, python-based, gpu-accelerated hyperspectral denoising package[END_REF], [START_REF] Rasti | Noise reduction in hyperspectral imagery: Overview and application[END_REF], [START_REF] Rasti | Image restoration for remote sensing: Overview and toolbox[END_REF] and noise reduction as a preprocessing step could improve the unmixing performance [START_REF] Rasti | How hyperspectral image unmixing and denoising can boost each other[END_REF]. We should note that unmixing techniques usually consider Gaussian noise and even can perform as a denoiser but are not as efficient as denoising techniques [START_REF] Rasti | How hyperspectral image unmixing and denoising can boost each other[END_REF]. Additionally, other types of noise, such as strip and sparse, can be removed by applying denoising before unmixing [START_REF] Rasti | Image restoration for remote sensing: Overview and toolbox[END_REF]. Alternatively, some methods were proposed for performing denoising and unmixing in a unified framework for boosting the performance of each other [START_REF] Ince | Simultaneous nonconvex denoising and unmixing for hyperspectral imaging[END_REF], [START_REF] Yang | Coupled sparse denoising and unmixing with low-rank constraint for hyperspectral image[END_REF]. As a result, we suggest applying a well-established noise reduction technique before estimating the abundance or characterizing the unmixing inverse problem using SNR-SD.

F. Unmixing Approaches w.r.t. Prior Knowledge of Endmembers

Unmixing techniques can be categorized in terms of prior knowledge of the endmembers. Depending on the prior knowledge available about endmembers, the unmixing problem can be divided into three main categories: (1) supervised unmixing, (2) unsupervised (blind) unmixing, and (3) Semi-supervised unmixing. In supervised unmixing, abundances are estimated by relying on known endmembers, whereas blind unmixing estimates both the endmembers and the abundances simultaneously. Semi-supervised unmixing relies on a library of endmembers that ideally contains the endmembers in the scene. We discuss each category separately. Fig. 7 demonstrates the graphical abstracts of these three linear unmixing categories. The taxonomy of the linear unmixing methods discussed in this paper is given in Fig. 9.

IV. SUPERVISED UNMIXING

In supervised unmixing, we assume that endmembers are known, and the abundance matrix needs to be estimated. In practice, the endmember can be measured in the field or laboratory. One can select them from a spectral library, however, due to the variations in the imaging setup, such a selection will not lead to a desirable abundance estimation. Alternatively, they can be extracted/estimated from data points. Geometrical approaches are often used at this step. Endmember extraction is not often easy since the captured data might not contain pure pixels to represent all the endmembers in the scene. Most endmember extraction methods rely on pure pixels or some pixels on the facet of the simplex. Geometrical approaches are less effective in no pure pixel scenarios. The processing chain of supervised unmixing is shown in Fig. 10. Note that we classify methods based on such sequential process (as opposed to simultaneous estimation of endmembers and abundances) as supervised unmixing since in the step of abundance estimation the endmembers are assumed to be known. The processing chain often includes three main steps; 1) Subspace projection: the data is projected into a subspace, 2) endmember extraction/estimation: Often a geometrical approach is used to extract the endmembers, and 3) Abundance estimation: abundances can be estimated using a least square (LS) or a neural network (NN) -based approach.

A. Subspace Projection and Estimation of number of Endmembers

Hyperspectral data often live in subspace having dimension much lower than the dimension of spectral bands defined by the sensor. Assuming r endmembers in the scene, the intrinsic/subspace dimension is r -1, i.e., the data points can be represented by r -1 linearly independent vectors or bases (in the case of orthogonal projections). Therefore, identifying such a subspace and projecting the data into that reduces the computational cost, memory consumption, and removes the noise and outliers.

The subspace projection is given by

Y = V T F + N, (22)
where columns V ∈ R p×r spans the subspace and columns of F are the projected spectral pixels. When V is semi orthogonal, i.e., V T V = I r , the project is orthogonal which is a common choice in HS subspace projection. In [START_REF] Heylen | Fully constrained least squares spectral unmixing by simplex projection[END_REF], both V and F should be estimated [START_REF] Rasti | Wavelet-Based Sparse Reduced-Rank Regression for Hyperspectral Image Restoration[END_REF], [START_REF] Rasti | Hyperspectral Feature Extraction Using Total Variation Component Analysis[END_REF]. Alternately, we can assume

[1] Y r = V T EA + V T N, (23)
Fig. 10. The processing chain of supervised unmixing. First, the data is projected into a subspace. Then, endmembers are extracted using a geometrical approach. If the endmembers are available, then the endmember extraction can be skipped. In the final step, abundances can be estimated using an ML/DLbased approach.

where Y r = V T Y. Reduction methods such as minimum noise fraction (MNF), principal component analysis (PCA), and noise adjusted principal component (NAPC) estimate V [START_REF] Rasti | Feature extraction for hyperspectral imagery: The evolution from shallow to deep (overview and toolbox)[END_REF]. Not that Y is usually not a square matrix; hence, singular value decomposition (SVD) can be used to estimate V [START_REF] Bioucas-Dias | Hyperspectral subspace identification[END_REF], [START_REF] Rasti | Hyperspectral subspace identification using SURE[END_REF]. Please see [START_REF] Ghamisi | Advances in hyperspectral image and signal processing: A comprehensive overview of the state of the art[END_REF], [START_REF] Rasti | Feature extraction for hyperspectral imagery: The evolution from shallow to deep (overview and toolbox)[END_REF] for more details on hyperspectral feature reduction.

Estimating the number of endmembers is not a trivial task. Unmixing approaches are vulnerable to this parameter and under or overestimating r may considerably affect the error of the models. We should note that sparse regression-based methods which do not use endmember bundles do not rely on the estimation of r. In the literature, this problem was addressed with alternative names such as hyperspectral subspace identification, intrinsic order selection, virtual dimension, and estimation of the number of spectrally distinct signal sources [START_REF] Bioucas-Dias | Hyperspectral subspace identification[END_REF]- [START_REF] Chang | Estimation of number of spectrally distinct signal sources in hyperspectral imagery[END_REF]. This problem can be addressed using eigenvaluebased detection techniques [START_REF] Chang | Estimation of number of spectrally distinct signal sources in hyperspectral imagery[END_REF]- [START_REF] Luo | Empirical automatic estimation of the number of endmembers in hyperspectral images[END_REF] or estimating the mean square errors [START_REF] Bioucas-Dias | Hyperspectral subspace identification[END_REF], [START_REF] Rasti | Hyperspectral subspace identification using SURE[END_REF]. Geometrical-based approaches were also proposed for endmember estimation [START_REF] Berman | Ice: a statistical approach to identifying endmembers in hyperspectral images[END_REF]- [START_REF] Ambikapathi | Hyperspectral data geometry-based estimation of number of endmembers using pnorm-based pure pixel identification algorithm[END_REF].

B. Endmember Extraction for Linear Unmixing

As can be seen in Fig. 11 (a), the endmembers are the vertices of the simplex enclosing the data points. Therefore, the geometrical concept inspires the mainstream of endmember extraction techniques. Geometrical approaches can be divided into two main groups; pure pixel-based approaches and approaches without pure pixel assumption. Fig. 11 compares the pure pixel with no pure pixel assumption. In Fig. 11 (b), the red circles demonstrate the positions of the missing pure pixels. In this case, the endmembers cannot be extracted and they should be estimated. This is feasible using geometrical approaches if enough pixels are located on the facets of the simplex. This becomes much more challenging for higher dimensional cases when noise and other nonlinearities are involved. Alternatively, a group of endmembers (endmember bundle) can be extracted based on the assumption of endmember variability caused by spectral variability as shown in Fig. 11 (c). Fig. 11 (d) shows the highly mixed scenario that geometrical approaches cannot cope with.

1) Pure Pixel Assumption: Many techniques consider pure pixel assumption for simplicity, while in real applications, pure pixels of some endmembers are often missing. The methods that rely on pure pixels for endmember extraction can be divided into three main groups: projections and extremes, simplex fitting methods, and multiple endmember extraction methods (endmember bundles). • Projections and Extremes: This group often searches for extremes by iteratively projecting data points. The vertices can be selected as the extreme points after iteratively projecting the data in a particular direction. For instance, Pixel Purity Index (PPI) [START_REF] Boardman | Mapping target signatures via partial unmixing of aviris data: in summaries[END_REF] scores the spectral vectors by projecting them onto a large set of random vectors (called skewers) and counting the number of times that each vector is an extreme point. Orthogonal subspace projection (OSP) [START_REF] Harsanyi | Hyperspectral image classification and dimensionality reduction: an orthogonal subspace projection approach[END_REF] and Vertex Component Analysis (VCA) [START_REF] Nascimento | Vertex component analysis: A˜fast algorithm to extract endmembers spectra from hyperspectral data[END_REF] selects endmembers iteratively by projecting the data into an orthogonal direction to the subspace spanned by the already selected endmember. • Simplex Fitting: The endmembers are assumed to be located at the vertices of the simplex enclosing the data points. Therefore, they can be extracted by maximizing the data simplex. N-FINDR [START_REF] Winter | N-FINDR: an algorithm for fast autonomous spectral end-member determination in hyperspectral data[END_REF] searches for pure pixels that form the largest simplex by gradually inflating a simplex inside the data. Simplex volume maximization (SiVM) [START_REF] Heylen | Fully constrained least squares spectral unmixing by simplex projection[END_REF] extracts the endmembers by iteratively maximizing the simplex volume using

arg max Er V (E r ) = arg max Er (-1) r • cmd (E r ) 2 r-1 (r -1)! , ( 24
)
where cmd is the Cayley-Menger determinant:

cmd

(E r ) = det        0 1 1 1 . . . 1 1 0 d 2 1,2 d 2 1,3 . . . d 2 1,r 1 d 2 2,1 0 d 2 2,3 . . . d 2 2,r 1 d 2 3,1 d 2 3,2 0 . . . d 2 3,r . . . . . . . . . . . . . . . . . . 1 d 2 r,1 d 2 r,2 d 2 r,3 . . . 0       
, and d i,j is the Euclidean distance between endmembers e i and e j . • Multiple Pure Pixels and Endmember Bundles: The multiple endmembers were also referred to as the endmember bundles. The endmember bundle and bundle unmixing were proposed in [START_REF] Bateson | A method for manual endmember selection and spectral unmixing[END_REF]- [START_REF] Bateson | Endmember bundles: a new approach to incorporating endmember variability into spectral mixture analysis[END_REF]. PCA was used to identify the (projected) endmembers using the extremes. More than one extreme was selected for materials. The bundle unmixing was formed using linear programming to determine the minimum, average, and maximum of the fractional abundances induced by endmember variability.

In [START_REF] Somers | Automated extraction of image-based endmember bundles for improved spectral unmixing[END_REF], endmember bundles were extracted from the observed data. Endmember variability is modeled using the endmember bundle, which will be discussed in more detail in semisupervised unmixing. Some techniques incorporate spatial information for extracting the endmembers [START_REF] Plaza | Spatial/spectral endmember extraction by multidimensional morphological operations[END_REF]- [START_REF] Shen | Spatial-spectral hyperspectral endmember extraction using a spatial energy prior constrained maximum simplex volume approach[END_REF]. These methods often utilize information from neighboring pixels in the endmember extraction. For instance, in [START_REF] Zortea | Spatial preprocessing for endmember extraction[END_REF], the similarity of the pixel and its neighbors was used to weigh the pixel and consequently adjust the simplex so that its vertices better represent the homogeneous regions in natural images. The primary assumption of these methods is that the neighboring pixels of pure pixels are also pure or have some purity level that can boost the endmember extraction. 2) No Pure Pixel Assumption: If there are no pure pixels, virtual endmembers can be successfully located at the vertices if enough data points are available on the facets of the data simplex. Here, the main idea is to minimize the simplex volume subject to a convex combination of the endmembers containing the data points. For instance, the minimum volume enclosing simplex (MVES) [START_REF] Chan | A convex analysis-based minimum-volume enclosing simplex algorithm for hyperspectral unmixing[END_REF], [START_REF] Craig | Minimum-volume transforms for remotely sensed data[END_REF] seeks the minimum simplex, which encloses the data points. Maximum Volume Inscribed Ellipsoid (MVIE) [START_REF] Lin | Maximum volume inscribed ellipsoid: A new simplex-structured matrix factorization framework via facet enumeration and convex optimization[END_REF] seeks the maximum volume ellipsoid contained within the convex hull of the data points. This volume is maximal when touching the facets of the data simplex. The contact points aid in estimating the virtual endmembers.

Alternatively, the minimum volume simplex analysis (MVSA) [START_REF] Li | Minimum volume simplex analysis: A fast algorithm to unmix hyperspectral data[END_REF] and the simplex identification via variable splitting and augmented Lagrangian (SISAL) [START_REF] Bioucas-Dias | A variable splitting augmented lagrangian approach to linear spectral unmixing[END_REF], [START_REF] Li | Minimum volume simplex analysis: A fast algorithm to unmix hyperspectral data[END_REF] form a non-convex minimization (or a non-concave maximization) problem which estimates endmembers by minimizing volume data simplex, subject to the ASC and ANC [START_REF] Bioucas-Dias | A variable splitting augmented lagrangian approach to linear spectral unmixing[END_REF], [START_REF] Li | Minimum volume simplex analysis: A fast algorithm to unmix hyperspectral data[END_REF]. Assuming

Q = E -1 r , SISAL solves Q = arg min Q -log|det(Q)| + λ1 T r hinge(QY r )1 n s.t. 1 T r Q = 1 T n Y -1 r . ( 25
)
where the hinge function (hinge(x) = max{-x, 0}) is the regularization allowing violations to the non-negativity constraint and makes SISAL robust to noise. SISAL and MVSA solve the same problem however, SISAL is more efficient since it uses an alternating direction method of multipliers (ADMM) [START_REF] Eckstein | On the douglas-rachford splitting method and the proximal point algorithm for maximal monotone operators[END_REF].

C. Abundance Estimation 1) Least Squares-based Approaches: In supervised unmixing, when the endmember E are known, the abundances should be estimated. This step is also known as inversion, referring to the inverse problem of estimating abundances when the endmembers are known. The presence of noise, inevitable errors in endmember estimation/extraction, and the physical constraints make such an inverse problem very challenging.

Unconstrained least squares unmixing (UCLSU) via the orthogonal subspace projection was proposed for abundance estimation [START_REF] Chang | Least squares subspace projection approach to mixed pixel classification for hyperspectral images[END_REF]. Non-negative constrained least squares unmixing (NCLSU) [START_REF] Settle | Linear mixing and the estimation of ground cover proportions[END_REF], [START_REF] Chang | Constrained subpixel target detection for remotely sensed imagery[END_REF] was proposed to estimate the abundances subjected to ANC. A weighted least square is suggested in [START_REF] Shimabukuro | The least-squares mixing models to generate fraction images derived from remote sensing multispectral data[END_REF] for estimating the abundances of multispectral remotely sensed data. There are several attempts to solve the least squares problem subjected to both ANC and ASC [START_REF] Shimabukuro | The least-squares mixing models to generate fraction images derived from remote sensing multispectral data[END_REF], [START_REF] Shimabukuro | Shade Images Derived from Linear Mixing Models of Multispectral Measurements of Forested Areas[END_REF]. In [START_REF] Shimabukuro | The least-squares mixing models to generate fraction images derived from remote sensing multispectral data[END_REF], quadratic programming was suggested for the constrained least squares. The first efficient algorithm was proposed in [START_REF] Heinz | Fully constrained least squares linear spectral mixture analysis method for material quantification in hyperspectral imagery[END_REF] and called fully constrained least squares unmixing (FCLSU). Later, fully constrained least squares unmixing by simplex projection was proposed [START_REF] Heylen | Fully constrained least squares spectral unmixing by simplex projection[END_REF] to derive a more efficient algorithm. A recursive algorithm minimizes the least square by performing orthogonal projections while holding the ASC and ANC. However, we should note that with the advances in graphical processing units (GPU), FCLSU can be efficiently solved using convex optimization techniques. Therefore, FCLSU is the most widely used method for abundance estimation given by:

 = arg min A 1 2 ||Y -EA|| 2 F s.t.A ≥ 0, 1 T r A = 1 T n . (26)
where the ||.|| F denotes the Frobenius norm. Problem ( 26) is convex and can be solved using any convex optimization solver. To promote spatial regularization, one can solve a penalized least square for abundance estimation given by

 = arg min A 1 2 ||Y -EA|| 2 F + λφ(A) s.t.A ≥ 0, 1 T r A = 1 T n , (27)
where the first term is the fidelity term. The penalty term is defined based on function φ applied to the unknown abundances. φ can be a total variation or sparsity-based function that usually incorporates spatial information. λ controls the trade-off between the fidelity and penalty terms. Problem [START_REF] Borsoi | A fast multiscale spatial regularization for sparse hyperspectral unmixing[END_REF] is convex (E is known) as long as the penalty term is convex and, therefore, it can be solved using any convex optimization, least squares, or quadratic programming solver. Problem [START_REF] Borsoi | A fast multiscale spatial regularization for sparse hyperspectral unmixing[END_REF] will lead to better estimation in terms of RMSE compared to [START_REF] Iordache | Collaborative sparse regression for hyperspectral unmixing[END_REF] for some values of λ. However, the selection of optimum λ is not trivial.

2) Shallow/Deep Neural Networks: In supervised unmixing networks, the endmembers or the decoder weights are fixed. Therefore, a regularized loss function can be given by L reg (y, ŷ) = L(y, ŷ) + λJ(a, W en ).

(28)

The regularizer J can be a combination of a Frobenius norm on W en with sparsity norms such as 1 norm on a. Autoencoders only utilize spectral information. Therefore, convolutional autoencoders were suggested for spectral unmixing to use both spatial and spectral information. Deep convolutional autoencoders were used for supervised hyperspectral unmixing in [START_REF] Khajehrayeni | Hyperspectral unmixing using deep convolutional autoencoders in a supervised scenario[END_REF], exploiting 3D convolutional filters.

In [START_REF] Rasti | UnDIP: Hyperspectral unmixing using deep image prior[END_REF], unmixing using deep image prior was proposed (UnDIP). Selection of the prior φ can be data dependent. Inspired by deep image prior [START_REF] Ulyanov | Deep image prior[END_REF], [START_REF] Ulyanov | Deep image prior[END_REF], UnDIP shows that the selection of the image prior φ in the minimization problem [START_REF] Borsoi | A fast multiscale spatial regularization for sparse hyperspectral unmixing[END_REF] can be shifted into optimizing a deep network's parameters. Then we have,

θ = arg min θ 1 2 ||Y -Ef θ (Z)|| 2 F s.t. Â = f θ (Z). ( 29
)
Z is the network input fixed throughout the learning. f θ is a deep network with parameters θ representing the weights and biases. UnDIP uses SiVM to extract E. Therefore, it is suitable for pure pixel-based scenarios. UnDIP uses an encoderdecoder network for f θ and (29) as the loss function and a softmax layer to enforce ASC and ANC. In [START_REF] Han | Deep halfsiamese networks for hyperspectral unmixing[END_REF], a two-stream network was proposed. The network contains an encoder in parallel with an autoencoder. First, bundle endmembers were extracted using VCA, and the corresponding abundances were estimated using the linear mixing model. Then, the encoder was trained based on the bundle endmembers and the abundances. The encoder shares the weights with the encoder of the autoencoder, whose decoder is nonlinear. A similar idea was also used in [START_REF] Hong | Endmember-guided unmixing network (egu-net): A general deep learning framework for self-supervised hyperspectral unmixing[END_REF].

V. UNSUPERVISED (BLIND) UNMIXING

In supervised unmixing, the processing chain is sequential. Usually, the estimation of abundance does not affect the estimation of endmembers due to the order in the processing chain. In unsupervised unmixing, we assume that both endmember and abundances are unknown and blind unmixing endmember and abundances are estimated simultaneously. We consider three major paradigms; 1) Least Squares-based Approaches 2) Shallow/Deep learning-based Approaches, and 3) Statistical-Based Approaches. Due to the inherent nonconvexity of blind unmixing methods, they are often vulnerable to initialization, and therefore they are always initialized using a geometrical endmember extraction approach.

A. Least Squares-based Approaches

Least squares-based approaches usually solve a minimization problem with respect to both E and A in the form of

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + λφ(A) + βψ(E) s.t.A ≥ 0, 1 T r A = 1 T n , 0 ≤ E ≤ 1 (30)
The main difference between the unmixing approaches in this group is the selection of the penalty functions φ and ψ. φ is often selected for incorporating spatial information and capturing the spatial correlation. The common choices for φ are TV and sparsity-promoting norms. ψ is usually a geometrical regularizer for incorporating geometrical information. 1) Selection of a Regularizer: The common choices for ψ are simplex volume terms and norms with simplex volume flavor, such as TV [START_REF] Berman | Ice: a statistical approach to identifying endmembers in hyperspectral images[END_REF]. Total variation penalty enforces the data simplex to have minimum volume [START_REF] Berman | Ice: a statistical approach to identifying endmembers in hyperspectral images[END_REF] ( Â, Ê) = arg min

A,E 1 2 ||Y -EA|| 2 F + βT V e (E) s.t.A ≥ 0, 1 T r A = 1 T n ( 31
)
where T V e (E) is a TV function given by

T V e (E) = r i,j=1 ||e i -e j || 2 2 = ||E(I r - 1 r 1 r 1 T r )|| 2 F . ( 32
)
See Appendix A for more details. We should mention that the method proposed in [START_REF] Berman | Ice: a statistical approach to identifying endmembers in hyperspectral images[END_REF] solves the problem after subspace projection using MNF. Minimum volume-constrained nonnegative matrix factorization (MVC-NMF) was proposed in [START_REF] Miao | Endmember extraction from highly mixed data using minimum volume constrained nonnegative matrix factorization[END_REF]. MVC-NMF suggests a minimum volume penalty for the constraint nonnegative least squares given by

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + βV ol 2 (E r ) s.t.A ≥ 0, 1 T r A = 1 T n , E ≥ 0. ( 33
)
Collaborative nonnegative matrix factorization (CoMNF) [START_REF] Li | Collaborative nonnegative matrix factorization for remotely sensed hyperspectral unmixing[END_REF] suggests a sparsity-promoting term on abundances and MV on the simplex given by

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + β E -m1 T r 2 F + λ r i=1 a (i) q 2 s.t. A ≥ 0, 1 T r A = 1 T n , (34)
where, 0 ≤ q ≤ 1, m contains the mean values of the spectral pixels, i.e., m = 1 n Y1 n , and a T (i) is the ith row of the matrix A. This term pulls the endmembers toward the center of mass. CoNMF uses both spatial and MV (geometrical) regularizers and solves the problem by projecting the data into a subspace. In [START_REF] Li | Robust collaborative nonnegative matrix factorization for hyperspectral unmixing[END_REF], Robust CoMNF (RCoNMF) was proposed, which utilizes a geometrical penalty that minimizes the distances between the endmembers to be estimated and the boundary pixels (P ∈ R p×r ). The main assumption is that the endmembers are close to the extremes of the data simplex (so-called boundary pixels). Hence, RCoNMF solves

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + β E -P 2 F + λ r i=1 a (i) 2 s.t. A ≥ 0, 1 T r A = 1 T n , (35)
where the boundary pixels, P, can be obtained using pure pixel-based endmember extraction. In [START_REF] Li | Robust collaborative nonnegative matrix factorization for hyperspectral unmixing[END_REF], VCA was used to estimate P. In [START_REF] Qian | Hyperspectral unmixing via l 1/2 sparsity-constrained nonnegative matrix factorization[END_REF], 0.5 -NMF was proposed which enforces q = r i=1 n j=1 |a ij | 0.5 sparsity penalty (q = 1/2) on the abundances. A multiplicative update rule [START_REF] Lee | Algorithms for non-negative matrix factorization[END_REF], [START_REF] Hoyer | Non-negative sparse coding[END_REF] was proposed to solve the 0.5 -NMF problem.

In such minimization problems, selecting the regulation parameters is not trivial. In [START_REF] Zhuang | Regularization parameter selection in minimum volume hyperspectral unmixing[END_REF], Nonnegative Matrix Factorization-Quadratic Minimum Volume (NMF-QMV) was proposed. NMF-QMV unifies the three geometrical regularizers mentioned above (i.e., TV, center, and boundary) in a single framework. NMF-QMV proposes the unified cost function using the MV-regularizers given by ( Â, Ê) = arg min

A,E 1 2 ||Y -EA|| 2 F + β EG -O 2 F s.t. A ≥ 0, 1 T r A = 1 T n , (36)
where the selection of G and O are given in Table II. We emphasize that those terms do not compute the volume of a simplex. They are convex penalties with volume flavor. A parameter selection technique was proposed in 36 by searching for an optimum β in a large interval. β optimum is when the estimated endmember simplex boundary is close to the data points boundary. Note that the spatial regularizer was ignored in [START_REF] Yu | Multi-stage convolutional autoencoder network for hyperspectral unmixing[END_REF], simplifying the selection technique. We should note that some blind unmixing techniques [START_REF] Zhuang | Regularization parameter selection in minimum volume hyperspectral unmixing[END_REF], [START_REF] Berman | Ice: a statistical approach to identifying endmembers in hyperspectral images[END_REF], [START_REF] Li | Collaborative nonnegative matrix factorization for remotely sensed hyperspectral unmixing[END_REF], [START_REF] Li | Robust collaborative nonnegative matrix factorization for hyperspectral unmixing[END_REF], solve the problems in a subspace. This has two main advantages: 1) reducing the computations, and 2) reducing noise and removing outliers. Hence, the endmembers are assumed to belong to an r -1 dimensional affine set, which best represents the data using an orthogonal projection. Note that the abundances do not change and can be estimated in a subspace. However, in real-world scenarios, the estimated endmembers in a subspace after reprojection into the original space may have negative values due to noise and nonlinearity since they cannot be bounded between zero and one in a subspace.

In [START_REF] Sigurdsson | Hyperspectral unmixing with lq regularization[END_REF], the spectral roughness penalty was proposed for imposing the smoothness to the endmembers, and the q norm was applied to the abundances to enforce sparsity on the abundances. The proposed cost function is given by

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + β 2 RE 2 F + λ n i=1 a i q + α 2 1 T r A -1 T n 2 2 s.t. A ≥ 0, (37)
where 0 ≤ q ≤ 1 and R ∈ R (p-1)×p is the first-order difference matrix given by

R =        -1 1 0 . . . 0 0 0 -1 1 0 • • • 0 . . . . . . . . . . . . . . . . . . 0 0 0 . . . -1 1        . ( 38
)
We should note that the roughness penalty does not incorporate geometrical information, and ASC constraint was relaxed and added to the cost function as a penalty. This makes the algorithm robust to the noise however, it adds one extra regularization parameter, α, which together with β and λ become very challenging to select, particularly in realworld applications. In [START_REF] Sigurdsson | Smooth spectral unmixing using total variation regularization and a first order roughness penalty[END_REF], the (2 dimensional) isotropic total variation was applied to abundances to promote piece-wise smoothness spatially. The cost function is given by

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + λ 1 r i=1 T V i (a (i) )+ λ 2 r i=1 a (i) q s.t. A ≥ 0, E ≥ 0, e (i) 2 2 = 1, i = 1, 2, ..., r, (39)
where 0 ≤ q ≤ 1 and T V i is the nonisotropic TV function given by

T V i (x) = (H h (x)) 2 + (H v (x)) 2 1 . (40)
H v = I ⊗ R and H h = R ⊗ I are the matrix operators to calculate the first order vertical and horizontal differences, respectively, for the vectorized image x. We should note that ASC was ignored in [START_REF] Rasti | Image restoration for remote sensing: Overview and toolbox[END_REF], and the endmember norm constraint (ENC) was applied to enforce a unit norm on the endmembers. ENC captures the spectral variability but has a big disadvantage, which does not let the algorithm distinguish the materials with similar endmembers with some scaling factors, such as tree and grass. In [START_REF] He | Total variation regularized reweighted sparse nonnegative matrix factorization for hyperspectral unmixing[END_REF], NMF subject to ASC was suggested using a reweighted 1 norm and the total variation penalty applied to the abundances. Many NMF-based blind unmixing approaches were proposed using the form of minimization [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF] [92]- [START_REF] Zhang | Spectral-spatial hyperspectral unmixing using nonnegative matrix factorization[END_REF].

A plug-and-play (PnP) priors framework was proposed in [START_REF] Zhao | A plug-and-play priors framework for hyperspectral unmixing[END_REF] for blind unmixing. Problem [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF], omitting the penalty and constraint on E, was solved using an ADMM. The penalized least squares step in the ADMM can be solved using a denoiser. Several denoisers, including a deep network, were suggested to solve that step. Additionally, the penalty term was considered not only applied to the abundances but also to the reconstructed data, i.e., φ(EA).

2) Archetypal Analysis Models: Several archetypal analysis (AA)-based [START_REF] Cutler | Archetypal analysis[END_REF] algorithms were proposed for blind unmixing. Let's assume that the endmembers are convex combinations of the spectral pixels, then using AA we have

( Â, B) = arg min A,B 1 2 ||Y -YBA|| 2 F s.t. B ≥ 0, 1 T n B = 1 T r , A ≥ 0, 1 T r A = 1 T n , (41)
where B ∈ R n×r and the columns of B belong to the simplex ∆ n . In [START_REF] Zouaoui | Entropic descent archetypal analysis for blind hyperspectral unmixing[END_REF], Entropic Descent Archetypal Analysis (EDAA) was proposed. EDAA uses an entropic descent algorithm to solve the AA problem [START_REF] Ince | Simultaneous nonconvex denoising and unmixing for hyperspectral imaging[END_REF]. A model selection technique suggested where the coherence µ = max i =j e T i e j and 1 residual ||Y -Ê Â|| 1 are minimized. This makes the solution robust to the initialization of B. The convex combination constraint on the endmembers is very strict. Therefore, in [START_REF] Zhao | Multiple endmembers based unmixing using archetypal analysis[END_REF], a variation of AA suggested which uses a relaxation factor δ and an r dimensional scaling vector s ∈ R r ( Â, B, Ŝ) = arg min

A,B,S 1 2 ||Y -YB(S I)A|| 2 F s.t. B ≥ 0, 1 T n B = 1 T r , A ≥ 0, 1 T r A = 1 T n , 1 -δ ≤ s ≤ 1 + δ (42)
where S = s1 T and denotes the Hadamard product. Additionally, Kernel AA (KAA) was suggested to take into account the endmember variability by multiple endmember extraction. The proposed KAA was computationally expensive, and therefore, in [START_REF] Zhao | Hyperspectral image unmixing based on fast kernel archetypal analysis[END_REF], the Nyström method is used to construct a low-rank approximation of the high-dimensional kernel matrix using K-means.

An 1 penalized variation of AA was proposed in [START_REF] Xu | L1 sparsity-constrained archetypal analysis algorithm for hyperspectral unmixing[END_REF] given by ( Â, B, Ŝ) = arg min

A,B,S 1 2 ||Y -PB(S I)A|| 2 F + ||A|| 1 s.t. B ≥ 0, 1 T n B = 1 T r , A ≥ 0, 1 T r A = 1 T n , 1 -δ ≤ s ≤ 1 + δ (43)
where the sparsity was enforced on the abundance and P ∈ R n×c , c >> r contains boundary pixels extracted by using PPI to cope with the endmember variability.

3) Extended Linear Mixing Models and Spectral variability: An Extended Linear Mixing Model (ELMM) was proposed in [START_REF] Veganzones | A new extended linear mixing model to address spectral variability[END_REF] to address the endmember variability given by

Y = ESA + N, A ≥ 0, , S ≥ 0, , E ≥ 0 ( 44
)
where S = s1 T and s contains the scaling factor of the endmembers. Nonnegative constraint least square (NCLS) was used to estimate the scaled abundances. Assuming the scaling factor is the same for all the endmembers and knowing that sum of the abundances is the scaling factor, therefore, the scaling factor was removed from abundances with a final normalization step. ELMM considers the spectral variation caused by illumination and topography. We should note that ELMM only considers uniform variations. In [START_REF] Drumetz | Blind hyperspectral unmixing using an extended linear mixing model to address spectral variability[END_REF], a pixeldependent ELMM was given by

y i = E i a i + n i , A ≥ 0, E i = E c i 1 T , a i ≥ 0, 1 T r a i = 1, (45)
where C ∈ R r×n contains the scaling factors in its column.

A minimization problem was proposed as

( Â, Ê, Ĉ) = arg min A,E,C 1 2 n i=1 ||y i -E i a i || 2 2 + λ 1 n i=1 ||E i -E c i 1 T || 2 F + λ 2 ψ(C) + λ 3 φ(A) (46)
s.t. A ≥ 0, 1 T r A = 1 T n
where E ∈ R p×r×n is a p × r × n tensor with E i located in band i and

ψ(C) = 1 2 (||RC|| 2 F + ||CR T || 2 F ) (47)
φ(A) = r j=1 ||H h (a (j) )|| 2 + ||H v (a (j) )|| 2 . ( 48
)
Generalized LMM (GLMM) was proposed in [START_REF] Imbiriba | Generalized linear mixing model accounting for endmember variability[END_REF] to address the nonuniform variation. GLMM is given by

y i = E i a i + n i , s.t. E i = E C i , a i ≥ 0, 1 T r a i = 1, (49)
where C i ∈ R p×r , i = 1, ..., n are scaling matrices dedicated to every pixel. The proposed minimization is given by

( Â, Ê, Ĉ) = arg min A,E,C 1 2 n i=1 ||y i -E i a i || 2 2 + λ 1 n i=1 ||E i -E C i || 2 F + λ 2 n i=1 ψ(C i ) + λ 3 φ(A) (50) s.t. A ≥ 0, 1 T r A = 1 T n ,
where C ∈ R p×r×n is p × r × n tensor with C i located in band i. In [START_REF] Borsoi | Improved hyperspectral unmixing with endmember variability parametrized using an interpolated scaling tensor[END_REF], a framework was suggested to estimate the scaling matrices C i using pure pixel endmember extraction techniques. A perturbed linear mixing model (PLMM) was proposed in [START_REF] Thouvenin | Hyperspectral unmixing with spectral variability using a perturbed linear mixing model[END_REF]. This simplifies the problem, however, the algorithm considerably relies on the presence of pure pixels and the extraction technique. The proposed model is given by

y i = E i a i +n i , s.t.E i = E+C i ≥ 0, E ≥ 0, a i ≥ 0, 1 T r a i = 1, (51)
where C i additive perturbation to model the spectral variability for every pixel. The proposed minimization is given by

( Â, Ê, Ĉ) = arg min A,E,C 1 2 n i=1 ||y i -E i a i || 2 2 + λ 1 2 ||AH d || 2 F + λ 2 ψ(E) + λ 3 n i=1 ||C i || 2 F s.t. A ≥ 0, 1 T r A = 1 T n , E ≥ 0, E i ≥ 0, i = 1, ..., n (52)
where ψ can be V ol 2 , T V e , or the Euclidean distance between E and some reference spectral signatures and H d is a matrix to calculate the deference of an abundance pixel with 4 nearest neighboring pixels (see [START_REF] Thouvenin | Hyperspectral unmixing with spectral variability using a perturbed linear mixing model[END_REF] for more details). PLMM was used in [START_REF] Thouvenin | Online unmixing of multitemporal hyperspectral images accounting for spectral vari-ability[END_REF], and GLMM was used in [START_REF] Borsoi | Kalman filtering and expectation maximization for multitemporal spectral unmixing[END_REF] for multitemporal unmixing.

B. Shallow/Deep Neural Networks

Shallow or deep neural networks can be used for blind unmixing, as shown in Fig. 8. The decoder is often shallow, and a linear layer is often used for linear unmixing. In blind unmixing networks, the decoder weights are the network's parameters representing the endmembers to be learned. Therefore, the regularized loss function is given by

L reg (y, ŷ) = L(y, ŷ) + λJ(a, W en , W de ) (53)
where J can be a combination of several functions, such as a Frobenius norm on W en with a sparsity norm such as 1 norm on a, and a Frobenius norm or a geometrical (MV induced) norm on W de . The loss function L might also be selected as a combination of several functions. For instance, EndNet [START_REF] Ozkan | Endnet: Sparse autoencoder network for endmember extraction and hyperspectral unmixing[END_REF] proposes a loss function with several terms, including a Kullback-Leibler divergence term, a SAD similarity, and a sparsity term, making the parameter selection very challenging. The blind unmixing networks can be divided into two main groups: 1) autoencoder-based networks and 2) algorithm unrolling-based networks.

1) Autoencoder-based networks: The majority of the blind unmixing networks are variations of autoencoder-based networks. The main components of such networks were explained in Section III-D. However, the encoder can be much deeper and more complicated. Five main types of AEs were considered for blind unmixing [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF]; 1) Sparse Nonnegative Autoencoders (SNAEs): SNAEs enforce sparsity on the bottleneck layer often via a regularization term in the loss, 2) Variational Autoencoders (VAEs): For VAEs, the decoder samples from a probability distribution model from the data and the encoded (hidden) features. The loss function uses a reconstructed loss plus a (KL-divergence) regularization term, ensuring that the encoded data stay close to a particular distribution (often a unit normal distribution), 3) Adversarial Autoencoders (AAEs): AAEs are based on generative adversarial networks (GANs). The network also uses a discriminator, and the AE's decoder often acts as a generator for the adversarial network. AAEs use a reconstructed plus an adversarial loss. 4) Denoising Autoencoders (DAEs): DAEs are trained based on partially corrupted input while minimizing the loss between the reconstructed and original (uncorrupted) inputs. Therefore, the network learns to remove the noise. 5) Convolutional Autoencoders (CAEs): CAEs use convolutional layers instead of fully connected layers. The weights are the convolutional filters that be learned throughout the training process. The convolutional filters capture the spatial correlation, which benefits the unmixing process. A comparative study was provided on the performance of AEs for hyperspectral unmixing in [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF]. The probabilistic generative model for spectral unmixing (PGMSU) was proposed in [START_REF] Shi | Probabilistic generative model for hyperspectral unmixing accounting for endmember variability[END_REF]. The network uses a VAE to generate the endmembers and a parallel encoder to generate the abundances.

SNSA [START_REF] Su | Stacked nonnegative sparse autoencoders for robust hyperspectral unmixing[END_REF] utilizes a stack of nonnegative sparse autoencoders. The last one performs the task of unmixing, and the others are exploited to improve the robustness concerning the outliers. Deep AutoEncoder Network (DAEN) [START_REF] Su | Daen: Deep autoencoder networks for hyperspectral unmixing[END_REF] exploits a stacked autoencoder to initialize a variational autoencoder which performs the unmixing task. An untied Denoising Autoencoder with Sparsity (uDAS) [START_REF] Qu | udas: An untied denoising autoencoder with sparsity for spectral unmixing[END_REF] exploits an additional denoising constraint on the decoder and a 2,1 sparsity promoting constraint on the decoder. Note that the autoencoder receives a spectral pixel at a time for training the network. Training a network based on a single spectrum at a time ignores the spatial information. The advantage of incorporating spatial information for spectral unmixing has been confirmed in the literature. Therefore, patch-wise or cube-wise CNN was proposed to utilize the spatial information. In [START_REF] Hua | Autoencoder network for hyperspectral unmixing with adaptive abundance smoothing[END_REF], an unmixing technique based on an autoencoder network is proposed that incorporates the spatial correlation between pixels by utilizing an adaptive abundance smoothing method. In [START_REF] Zhang | Hyperspectral unmixing via deep convolutional neural networks[END_REF], patch-wise or cube-wise convolutional autoencoders were proposed to incorporate the spatial information. In [START_REF] Palsson | Spectral-spatial hyperspectral unmixing using multitask learning[END_REF], parallel autoencoders were applied to spectral patches to exploit the spatial information. A CAE was proposed in [START_REF] Palsson | Convolutional autoencoder for spectral-spatial hyperspectral unmixing[END_REF]. The proposed network uses two convolutional layers as the decoder and a convolutional layer with linear activation as the decoder. ASC and ANC were enforced using a softmax, and SAD was used as the loss function. The cycle-consistency unmixing network [START_REF] Gao | Cycunet: Cycle-consistency unmixing network by learning cascaded autoencoders[END_REF] utilizes two convolutional autoencoders, which are cascaded and performed cyclically. The proposed loss function contains two terms for spectral reconstruction and one for abundance reconstruction to incorporate highlevel semantic information. In [START_REF] Ghosh | Hyperspectral unmixing using transformer network[END_REF], a combination of a CAE and a transformer was proposed for unmixing. The attention module of the transformer captures the nonlocal contextual feature dependencies between the image patches. In [START_REF] Qiwen | Adversarial autoencoder network for hyperspectral unmixing[END_REF], VCA and FCLSU were first used to drive prior abundances from superpixels. Then, the adversarial training guides the encoder of an AAE to transfer the spatial information into the network by matching the abundances. A VAE was used in [START_REF] Borsoi | Deep generative endmember modeling: An application to unsupervised spectral unmixing[END_REF] to address the spectral variability. First, a deep generative endmember model was learned using the pure pixel extracted from the data. Then, an optimization was suggested and solved using the endmember model as a prior. A deep generative model was proposed in [START_REF] Shi | Deep generative model for spatial-spectral unmixing with multiple endmember priors[END_REF], which uses three types of encoders, i.e., CNN, convolutional graph networks, and selfattention-based networks. A discriminator was used to enforce the generated endmembers to have a distribution similar to the endmember models either extracted from the data or selected from a library.

In AEs, a deeper encoder does not necessarily improve the performance of AE-based unmixing, however, it helps to enlarge the receptive field [START_REF] Palsson | Blind hyperspectral unmixing using autoencoders: A critical comparison[END_REF]. Therefore, in [START_REF] Yu | Multi-stage convolutional autoencoder network for hyperspectral unmixing[END_REF], a Multi-Stage CAE Network (MSNet) was suggested, which uses a multi-stage framework. Every stage uses a shallow encoderdecoder architecture. The input data are spatially downsampled at every stage. The features at the next stage were concatenated with the previous stage's input to preserve the abundance consistency of different stages.

A drawback of the many DL-based approaches is the absence of geometrical information. Most of them rely on SAD and/or MSE for training. As we already discussed and it has been proven [START_REF] Lin | Identifiability of the simplex volume minimization criterion for blind hyperspectral unmixing: The no-pure-pixel case[END_REF], geometrical information could be very beneficial without pure pixels. To incorporate geometrical and spatial information, in [START_REF] Rasti | Misicnet: Minimum simplex convolutional network for deep hyperspectral unmixing[END_REF], a minimum simplex convolutional network (MiSiCNet) was proposed. In [START_REF] Rasti | Misicnet: Minimum simplex convolutional network for deep hyperspectral unmixing[END_REF], it was shown that the minimization problem

( Â, Ê) = arg min A,E 1 2 ||Y -EA|| 2 F + β E -m1 T r 2 F + λφ(A) s.t. A ≥ 0, 1 T r A = 1 T n , 0 ≤ E ≤ 1 ( 54
)
can be turned into

( θ1 , θ2 ) = arg min θ1,θ2 1 2 ||Y -θ 2 f θ1 (Z)|| 2 F + λ θ 2 -m1 T r 2 F s.t. Ê Â = θ2 f θ1 (Z) (55)
where Ŷ = Ê Â, θ 1 and θ 2 are the parameters of the encoder and the decoder of the network, respectively. Therefore, the minimization problem ( 54) can be shifted to the minimization of the parameters of a deep network f with a loss function given by:

L(Y, Ŷ, Ê, m) = 1 2 ||Y -Ŷ|| 2 F + λ Ê -m1 T r 2 F , ( 56
)
where f is an encoder-decoder network and Ê is the weight of the decoder. MiSiCNet uses a softmax layer to enforce both the ASC and the ANC.

2) Algorithm Unrolling-based Networks: The main idea behind algorithm unrolling-based networks is to incorporate model-based with learning-based methods. In this group, an iterative algorithm (such as ADMM) is considered to solve a mixture model. Then, the iterations of the algorithm will be unrolled and implemented using neural networks. The updating equations will be substituted using NNs. Sequential repetitions of such networks will be trained using a loss function such as MSE. An advantage of such networks is that the tuning parameters can be trainable. SUnSAL with ADMM implementation was unrolled in [START_REF] Zhou | Admm-based hyperspectral unmixing networks for abundance and endmember estimation[END_REF] with an extra linear layer augmented to the network to implement blind unmixing (ADMMNet). The sequential networks in [START_REF] Zhou | Admm-based hyperspectral unmixing networks for abundance and endmember estimation[END_REF] were tested with and without parameter sharing. In [START_REF] Qian | Spectral mixture model inspired network architectures for hyperspectral unmixing[END_REF], 1 -NMF was considered for blind unmixing. The iterative softthresholding algorithm [START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF], was unrolled to estimate the abundances. An extra linear layer was used for the endmember estimation and added to the sequential networks. The scenarios of the networks with and without sharing parameters were also considered in [START_REF] Qian | Spectral mixture model inspired network architectures for hyperspectral unmixing[END_REF]. A solution to the p -NMF was unrolled for blind unmixing in [START_REF] Xiong | Snmf-net: Learning a deep alternating neural network for hyperspectral unmixing[END_REF]. Instead of using a linear layer, the endmember estimation step was also unrolled. Therefore, the sequential network was constructed of endmember networks followed by abundance networks. For the abundance and endmember networks, the generalized shrinkage thresholding solution given in [START_REF] Zuo | A generalized iterated shrinkage algorithm for non-convex sparse coding[END_REF] and the projected gradient descent were unrolled, respectively.

Finally, we should note that often blind NN-based unmixing can be turned into supervised unmixing by fixing the weights of the decoders using extracted endmembers.

C. Statistical-based Approaches

Most recent advanced techniques developed for blind unmixing fall into the first two categories. There are, however, other statistical-based approaches that cannot be reformulated as the constraint penalized least squares. The Bayesian approach is widely used for blind unmixing. The maximum a posteriori (MAP) estimator and the minimum mean squared error (MMSE) estimator were used in [START_REF] Parra | Unmixing hyperspectral data[END_REF] and [START_REF] Dobigeon | Joint bayesian endmember extraction and linear unmixing for hyperspectral imagery[END_REF], respectively. The Bayesian approaches using the MAP estimator can be reformulated as the penalized least square in the form of [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF] when the random noise is assumed to be Gaussian and uncorrelated, i.e., the noise covariance matrix is diagonal [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], [START_REF] Elad | Analysis versus synthesis in signal priors[END_REF]. This assumption simplifies the inverse problem and its solution. Independent Component Analysis (ICA) [START_REF] Jutten | Blind separation of sources, part 1: An adaptive algorithm based on neuromimetic architecture[END_REF], as a solution to blind source separation, was proposed for blind unmixing [START_REF] Bayliss | Analyzing hyperspectral data with independent component analysis[END_REF], [START_REF] Chen | Independent component analysis for remote sensing study[END_REF]. The abundances are dependent (ASC), therefore, ICA is not a suitable solution to blind unmixing [START_REF] Nascimento | Does independent component analysis play a role in unmixing hyperspectral data?[END_REF]. Dependent component analysis (DECA) was proposed for blind unmixing of highly mixed datasets [START_REF] Nascimento | Hyperspectral unmixing algorithm via dependent component analysis[END_REF], [START_REF] Nascimento | Hyperspectral unmixing based on mixtures of dirichlet components[END_REF]. DECA assumes a mixture of Dirichlet densities as the prior for abundances. Dirichlet density satisfies ANC and ASC. The Bayesian estimation problems are often hard to solve, and the estimation of unknown parameters is computationally expensive. In [START_REF] Dobigeon | Bayesian separation of spectral sources under non-negativity and full additivity constraints[END_REF], a Gibbs sampling strategy was proposed to decrease the complexity of the proposed hierarchical Bayesian model. A similar Bayesian approach was proposed in [START_REF] Dobigeon | Joint bayesian endmember extraction and linear unmixing for hyperspectral imagery[END_REF] for spectral mixture analysis (BayesianSMA) from Raman spectroscopy.

VI. SEMISUPERVISED (LIBRARY-BASED) UNMIXING

When neither the pure pixels nor the sufficient spectra on the facets of the data simplex are available, the endmembers can not be successfully extracted/estimated, which results in poor abundance estimations. Alternatively, blind unmixing techniques can be used. However, blind unmixing is a nonconvex problem; therefore, a good initialization of endmembers using endmember extraction techniques will often benefit in finding a better solution. Additionally, in highly mixed scenarios, blind unmixing methods often fail due to the large sets of solutions that fit the data. Fig. 11 (d) demonstrates that in highly mixed scenarios, finding the real simplex is not trivial. Therefore, semi-supervised unmixing was suggested, relying on an endmember library. Here, we should note that we categorize all techniques that entirely or partially rely on a spectral library as semi-supervised unmixing techniques. Therefore, before discussing the semisupervised techniques, we discuss the library selection and construction.

Library selection or construction is a crucial step to the success of a semisupervised unmixing method. Blind selection of a library without extra attention and some processing steps will lead to poor results for semisupervised unmixing. There are two major paradigms in semisupervised unmixing i) Multiple Endmember Spectral and Mixture Analysis and ii) Sparse unmixing. The former was designed to address the spectral variability using the endmember variability. Therefore, the spectral library is designed to represent the variability of the endmembers. The latter seeks a sparse solution relying on the library, and therefore, the high correlation of the library endmembers avoids the sparse solution. In both paradigms, the library must well represent the materials in a sense, i.e., it must contain all the endmembers of the materials in the scene. Generally, a library can be obtained using either of the following ways [START_REF] Borsoi | Spectral variability in hyperspectral data unmixing: A comprehensive review[END_REF], [START_REF] Somers | Chapter 17endmember library approaches to resolve spectral mixing problems in remotely sensed data: Potential, challenges, and applications[END_REF].

1) In situ field or/and laboratory measurements: endmember libraries can be built using field or/and laboratory measurements [START_REF] Roberts | Mapping chaparral in the santa monica mountains using multiple endmember spectral mixture models[END_REF]. Creating such libraries has several major drawbacks. It is hard, time-consuming, expensive, and might be sensor-dependent. Due to the different measuring conditions, systems, or instruments, there are often mismatches and scaling discrepancies between the endmembers from the library and the observed data [START_REF] Tompkins | Optimization of endmembers for spectral mixture analysis[END_REF], [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF]. A list of available spectral libraries was given on https://specchio.ch/. 2) Construction using observed data: Extracting multiple pure pixels from the observed datasets representing every material in the observed dataset. These spectra can be further clustered to represent the endmember bundles. These libraries are data-dependent and often fail in highly mixed scenarios since there are not enough pure pixels in the dataset [START_REF] Somers | Automated extraction of image-based endmember bundles for improved spectral unmixing[END_REF], [START_REF] Quintano | Multiple endmember spectral mixture analysis (mesma) to map burn severity levels from landsat images in mediterranean countries[END_REF], [START_REF] Drumetz | Hyperspectral image unmixing with endmember bundles and group sparsity inducing mixed norms[END_REF]. 3) Construction using physical models: Radiative transfer models can be used to generate endmember libraries [START_REF] Tits | The potential and limitations of a clustering approach for the improved efficiency of multiple endmember spectral mixture analysis in plant production system monitoring[END_REF]. A radiative transfer model is defined as a function of physicochemical parameters. Therefore, different instances of endmembers can be obtained by varying physicochemical parameters. For instance, PROSPECT model [START_REF] Jacquemoud | Prospect: A model of leaf optical properties spectra[END_REF] for vegetation and Hapke model [START_REF] Hapke | Bidirectional reflectance spectroscopy: 1. theory[END_REF] for reflectance modeling of intimate mixtures and densely packed grains or particles can be used. These libraries are sensor-and data-independent. However, they depend on the availability and accuracy of a model for the target materials. Additionally, selecting representative samples from such complex models is a big challenge. A combination of the libraries mentioned above also may be used. Alternatively, generative machine learning models can be used to augment a library. In [START_REF] Borsoi | Deep generative models for library augmentation in multiple endmember spectral mixture analysis[END_REF], a variational autoencoder was used for library augmentation.

A. Multiple Endmember Spectral and Mixture Analysis

As we discussed, endmember variability is a challenge in hyperspectral unmixing. A solution is to use a library containing endmembers' variations for materials. Then, the unmixing problem can be formulated as a semisupervised problem relying on such a dictionary. The pioneer method was proposed in [START_REF] Roberts | Mapping chaparral in the santa monica mountains using multiple endmember spectral mixture models[END_REF] called multiple endmember spectral and mixture analysis (MESMA). MESMA was proposed to address the endmember variability. Assuming a structured library, D = [D 1 , D 2 , ..., D r ] contains the endmember bundles, D i 's, for all the materials, MESMA allows different and scaled endmembers for every pixel solving the problem âi , Êi = arg min ai,Ei

1 2 ||y i -E i a i || 2 2 s.t. E i ∈ D, a i ≥ 0, 1 T r a i = 1. (57)
This is a combinatorial problem, and MESMA searches for all combinations of endmembers (from every endmember bundle) and selects the one with the lowest reconstruction error. Endmember candidates are sometimes called models. This is solved for all the pixels. MESMA is a highly computationally demanding algorithm without any stopping criteria Π i m i (m i 's are dimensions of endmember bundles) times FCLSU must be computed for every pixel. However, this can be performed in parallel for all the pixels, it is computationally costly, and therefore, several algorithms have been proposed to reduce the computational burden. As a result, instead of FCLSU, the sum to one constrained least square unmixing (SCLSU) was used with a threshold for the error (RE). Note that ignoring the nonnegativity constraint leads to SCLSU, which has a closed-form solution and can be implemented efficiently. Then, solutions that are positive and less than one were selected [START_REF] Roberts | Mapping chaparral in the santa monica mountains using multiple endmember spectral mixture models[END_REF], [START_REF] Heylen | Hyperspectral unmixing with endmember variability via alternating angle minimization[END_REF]. Additionally, the number of iterations and the threshold value can be changed to decrease the computations, but it can highly affect the performance [START_REF] Philip | A comparison of error metrics and constraints for multiple endmember spectral mixture analysis and spectral angle mapper[END_REF].

In [START_REF] Combe | Analysis of omega/mars express data hyperspectral data using a multiple-endmember linear spectral unmixing model (melsum): Methodology and first results[END_REF], Unconstrained Least Squares (UCLS) were used to speed up the search for the solutions to the combinatorial problem.

Alternatively, the size of the library can be decreased [START_REF] Roberts | Evaluation of the potential of hyperion for fire danger assessment by comparison to the airborne visible/infrared imaging spectrometer[END_REF]- [START_REF] Philip | Endmember selection for multiple endmember spectral mixture analysis using endmember average rmse[END_REF]. In [START_REF] Heylen | Hyperspectral unmixing with endmember variability via alternating angle minimization[END_REF], the angle between a pixel and the projections into hyperplanes of the fixed point (selected endmembers) and the endmember candidates were used to find the optimal endmember. The minimum angle is equivalent to the minimum error, and therefore, the SCLSU problem (note that projection into an affine set or a hyperplane only satisfies ASC) was shifted to more efficient angle computations. FCLSU was used to estimate the abundances using the selected endmembers (the optimal model), and the result with minimum reconstruction error was selected. On the other hand, in [START_REF] Borsoi | Deep generative models for library augmentation in multiple endmember spectral mixture analysis[END_REF], a deep generative model using variational autoencoders was suggested for augmenting the library endmembers to address the library mismatch. First, VAE was used to augment the library, and then MESMA was used for unmixing. However, this makes the algorithm computationally very expensive. A similar but more efficient approach was proposed in [START_REF] Tits | A geometric unmixing concept for the selection of optimal binary endmember combinations[END_REF], which only uses one unmixing step.

B. Sparse Unmixing

The first efficient idea was based on sparse and redundant modeling [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF] and sparse regression, and therefore, it is often referred to as sparse unmixing. We should note that semisupervised unmixing refers to all the unmixing methods relying on a library and is not limited to sparse unmixing methods.

In sparse unmixing, the fractional abundances are estimated using sparse regression techniques. These methods describe each spectrum as a sparse linear combination of the elements of a rich library of pure spectra.

The main idea of sparse unmixing was first suggested in [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF]. If a well-designed dictionary is available, the sparse regression can be used to estimate the abundances given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ||X|| 1 s.t. X ≥ 0, 1 T m X = 1 T n . (58)
An ADMM-based algorithm was proposed to solve the problem ( 58), and therefore it was called sparse unmixing by variable splitting and augmented Lagrangian (SUnSAL). In [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF], also the constraint SUnSAL (CSUnSAL) proposed given by X = arg min

X ||X|| 1 s.t. ||Y -DX|| F < δ, X ≥ 0, 1 T m X = 1 T n , (59)
and was solved using an ADMM-based algorithm. However, it is suggested to use SUnSAL without ASC due to the conflict with 1 [START_REF] Sigurdsson | Hyperspectral unmixing with lq regularization[END_REF]. Additionally, ASC was found to be a rigorous constraint that often does not occur in the real world due to noise and signature variability [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF]. Therefore, SUnSAL often refers to the problem (58) without ASC, and problem [START_REF] Bajjouk | Quantification of subpixel cover fractions using principal component analysis and a linear programming method: Application to the coastal zone of roscoff (france)[END_REF] refers to FCLS version of SUnSAL which can be seen in the arXiv version of SUnSAL [START_REF] Bioucas-Dias | Alternating direction algorithms for constrained sparse regression: Application to hyperspectral unmixing[END_REF]. In [START_REF] Nascimento | Parallel hyperspectral unmixing on gpus[END_REF], problem ( 58) is suggested for implementation on a GPU when λ = 0. In [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF], the nonconvex 0 sparse regression X = arg min

X ||X|| 0 s.t. ||Y -DX|| F < δ, X ≥ 0, (60)
were compared with SUnSAL and CSUnSAL. The nonnegative 0 sparse regression problem (60) was solved using a greedy algorithm called orthogonal matching pursuit [START_REF] Pati | Orthogonal matching pursuit: Recursive function approximation with applications to wavelet decomposition[END_REF]. It was shown that SUnSAL outperforms the other two sparse regressions. Additionally, the experiments in [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF] confirmed the advantage of ANC. The success of sparse unmixing depends on the following assumptions; i) the library must well represent the materials in a sense, i.e., it must contain all the endmembers of the materials in the scene; ii) sparse solutions of the underdetermined linear systems exist. The latter depends on the sparseness degree of the observed spectra over the library (i.e., the number of library members representing the spectral pixels) and the coherence of library members. The coherence of a library can be measured using mutual coherence [START_REF] Donoho | Optimally sparse representation in general (nonorthogonal) dictionaries via 1 minimization[END_REF] (or restricted isometric constants (RICs) [START_REF] Candes | Stable signal recovery from incomplete and inaccurate measurements[END_REF]). A higher mutual coherence requires a lower sparseness degree to achieve a sparse solution. In sparse unmixing, the mutual coherence of the spectral library is very high [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF], which is the main drawback, however, the sparseness degree is low (i.e., pixels are mixtures of a few pure spectra). To decrease the mutual coherence, the library must be pruned before applying the sparse regression or throughout the estimation process. Alternatively, the sparsityinducing prior can be better defined by incorporating spatial information and/or applying segmentation, localization, and group sparsity, which will be discussed in detail in this section.

1) Selection of a Regularizer: In sparse unmixing, many attempts have been made to improve the sparse representation of the abundances over a spectral library by selecting different regularizers. In [START_REF] Iordache | Total variation spatial regularization for sparse hyperspectral unmixing[END_REF], SunSAL was improved by incorporating spatial information and applying a TV penalty on abundances (SUnSAL-TV) given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ 1 ||X|| 1 + λ 2 T V ni (X) s.t. X ≥ 0. (61)
where T V ni (X) is a nonisotropic TV function given by

T V ni(X) = ||XH T h || 1 + ||XH T v || 1 (62)
One of the main challenges of sparse regression is the high mutual coherence of the spectral library, which might lead to poor abundance estimation [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF]. Therefore, in [START_REF] Iordache | Hyperspectral unmixing with sparse group lasso[END_REF], group SUNSAL (GSUNSAL) was suggested by dividing the library into L groups (structured dictionary) and adding a group LASSO (Least Absolute Shrinkage and Selection Operator) [START_REF] Meier | The group lasso for logistic regression[END_REF] X = arg min

X 1 2 ||Y -DX|| 2 F + n i=1 L g=1 λ g 1 ||x g i || 2 + λ 2 ||X|| 1 s.t. X ≥ 0. (63)
where x g is the gth subvector of abundances. As we mentioned, 1 norm conflicts with ASC and therefore q norm also investigated for sparse unmixing [START_REF] Chen | Sparse hyperspectral unmixing based on constrained lp -l2 optimization[END_REF] X = arg min

X 1 2 ||Y -DX|| 2 F + λ n i=1 ||x i || q q s.t. X ≥ 0, 1 T m X = 1 T n , 0 < q < 1 (64)
where

||x|| q q m i=1 |x i | q , (65)
Collaborative sparse unmixing [START_REF] Iordache | Collaborative sparse regression for hyperspectral unmixing[END_REF] applies 2,1 (i.e., the sum of 2 on the abundances) to promote the sparsity on the abundances.

X = arg min X 1 2 ||Y -DX|| 2 F + λ m i=1 ||x (i) || 2 s.t. X ≥ 0. (66)
Collaborative sparse unmixing using 0 norm was proposed in [START_REF] Shi | Collaborative sparse hyperspectral unmixing using l 0 norm[END_REF] given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ||X|| 2,0 s.t. X ≥ 0. (67)
where

||X|| 2,0 = m i=1 Id(||x (i) || 2 > 0)
and Id is the Heaviside step indicator function, i.e., it is one when the argument is positive. In [START_REF] Drumetz | Hyperspectral image unmixing with endmember bundles and group sparsity inducing mixed norms[END_REF], the endmember variability was considered using a dictionary of endmember bundles generated from the data. Unlike the methods mentioned above, in [START_REF] Drumetz | Hyperspectral image unmixing with endmember bundles and group sparsity inducing mixed norms[END_REF], the dictionary is extracted from data points. In [START_REF] Drumetz | Hyperspectral image unmixing with endmember bundles and group sparsity inducing mixed norms[END_REF], a framework of sparse regression was proposed given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ||X|| g,q,d s.t. X ≥ 0, 1 T m X = 1 T n , (68)
where the grouped sparse norm is defined as

||X|| g,q,d n i=1 ( L g=1 ||x g i || d q ) 1 d . (69)
Three cases are compared, i.e., {q = 1, d = 2}, {q = 1, 0 < d < 1}, and {q = 2, d = 1}. We should note that ASC was used in the endmember bundles method [START_REF] Drumetz | Hyperspectral image unmixing with endmember bundles and group sparsity inducing mixed norms[END_REF], unlike the structured dictionary methods [START_REF] Iordache | Collaborative sparse regression for hyperspectral unmixing[END_REF], [START_REF] Iordache | Hyperspectral unmixing with sparse group lasso[END_REF].

In [START_REF] Wang | Hyperspectral unmixing using double reweighted sparse regression and total variation[END_REF], a spatial-spectral weighted 1 norm used for (61) given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ 1 ||W SX|| 1 + λ 2 T V ni (X) s.t. X ≥ 0. ( 70
)
where S = s1 T , s i = (||x (i) || 1 + ) -1 , and W = (X+ ) -1 are the weight matrices updated in every iterations. The weights further promote the sparsity on X.

Local collaborative sparse unmixing (LCSU) [START_REF] Zhang | Hyperspectral unmixing based on local collaborative sparse regression[END_REF] and a spectral-spatial weighted sparse unmixing (S 2 WSU) [START_REF] Zhang | Spectral-spatial weighted sparse regression for hyperspectral image unmixing[END_REF] further incorporate spatial information by adapting the collaborative sparse unmixing and the reweighted 1 norm locally to neighboring pixels, respectively. Indeed, S 2 WSU uses problem (70) locally and when λ 2 = 0.

In [START_REF] Giampouras | Simultaneously sparse and low-rank abundance matrix estimation for hyperspectral image unmixing[END_REF], the weighted nuclear norm was added to the weighted 1 norm to utilize the low-rankness of the abundances. The proposed sparse regression is given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ 1 ||W X|| 1 + λ 2 ||X|| t, * s.t. X ≥ 0. ( 71
)
where the weighted nuclear norm is given by

||X|| t, * = r i=1 t i σ i (X), (72)
and σ i (X) are the singular values of X and the weights are given by t i = (σ i (X) + ) -1 . In [START_REF] Li | Hyperspectral sparse unmixing with spectral-spatial low-rank constraint[END_REF], the spatial-spectral weighted 1 was used with the weighted nuclear norm. In [START_REF] Han | Hyperspectral unmixing via nonconvex sparse and lowrank constraint[END_REF], the weighted Schatten p-norm [START_REF] Xie | Weighted schatten p -norm minimization for image denoising and background subtraction[END_REF] and a collaborative norm was suggested to induce the low rankness and the sparsity, respectively, given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ 1 m i=1 ||x (i) || q 2,q + λ 2 ||X|| q t,Sp s.t. X ≥ 0, 0 < q ≤ 1. (73)
where

||x (i) || q 2,q m i=1 ||x (i) || q 2 , (74)
||X|| q t,Sp r i=1 t i σ d i (X). (75)
2) Segmentation/Super-pixel and Spectral Variability: Some methods apply segmentation prior to sparse unmixing for grouping. In this way, they remove spectral variation and capture spatial information. A fast Multiscale Sparse Unmixing Algorithm (MUA) was proposed in [START_REF] Borsoi | A fast multiscale spatial regularization for sparse hyperspectral unmixing[END_REF]. MUA first applies segmentation techniques such as a binary partition tree (BPT), simple linear iterative clustering (SLIC), or K-means to group pixels. Then SUnSAL was utilized for unmixing the average pixels of every segment. The coarse fractional abundance matrix estimated using SUNSAL (X seg. ) is used as a prior for the final sparse regression given by

X = arg min X 1 2 ||Y -DX|| 2 F + λ 1 ||X|| 1 + λ 2 ||X -X seg. || 2 2 s.t. X ≥ 0. ( 76
)
This idea is also called "super-pixel" in the literature, and numerous methods have been proposed using this terminology. In [START_REF] Ince | Superpixel-based graph laplacian regularization for sparse hyperspectral unmixing[END_REF], SLIC was utilized for the pixel grouping and pixelbased sparse unmixing was performed using superpixel-based graph Laplacian regularization. In [START_REF] Zhang | Spectral variability augmented sparse unmixing of hyperspectral images[END_REF], the spectral variability was considered using a spectral variability dictionary. The dictionary was obtained by applying PCA on an endmember library extracted from the dataset using a pure pixel-based endmember extraction [START_REF] Mei | Spatial purity based endmember extraction for spectral mixture analysis[END_REF]. The PCs with high variance were used to create the spectral variability library.

3) Library pruning: Library pruning is a necessary preprocessing step in sparse unmixing to reduce the mutual coherence of the spectral library. The pruning step considerably improves the computational time. However, it risks losing endmembers that are scaled versions of each other. Nevertheless, in real-world applications, the library pruning must be done cautiously, particularly when the SNR is low, and the number of endmembers is high [START_REF] Iordache | Potential and limitations of band selection and library pruning in sparse hyperspectral unmixing[END_REF]. A simple pruning often used before applying sparse unmixing is to compute SAD. Then, spectra with small angles are removed from the library to reduce the mutual coherence. A pruning strategy was proposed in [START_REF] Iordache | Dictionary pruning in sparse unmixing of hyperspectral data[END_REF]. The library was projected into the data subspace, and spectra with high projection errors (the normalized Euclidean distance) were removed from the library. Arguably, data-dependent library pruning is more suitable for sparse unmixing. In [START_REF] Iordache | Musiccsr: Hyperspectral unmixing via multiple signal classification and collaborative sparse regression[END_REF], a library pruning-based sparse unmixing called multiple signal classification collaborative sparse unmixing (MUSIC-CSR) was proposed. The library was pruned (before applying the collaborative sparse unmixing) using an orthogonal projection (P ⊥ = I -VV T ) [START_REF] Schmidt | Multiple emitter location and signal parameter estimation[END_REF] where V is the bases of the Y's subspace obtained by HySime [START_REF] Bioucas-Dias | Hyperspectral subspace identification[END_REF] to reduce the noise effect. The normalized projection error of each library endmember ||P ⊥ d i || 2 /||d i || 2 was sorted by increasing order. Hence, the first r defines the reduced library. A similar idea was used in [START_REF] Xu | Simultaneously multiobjective sparse unmixing and library pruning for hyperspectral imagery[END_REF], which formulates the library reduction into a multiobjective optimization w.r.t. a binary vector representing the support of the abundances matrix.

4) Library Mismatch: As mentioned, the endmember library is fixed in the problems and methods discussed, and the ultimate goal is abundance estimation. However, even a pruned and well-selected spectral library cannot represent all the endmembers of materials in a real-world dataset. There are several factors, such as noise, atmospheric effects, illumination variations, and the intrinsic variation of materials, which may change the endmembers and induce scaling factors into the endmembers in the scene compared to the ones from the library [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF], [START_REF] Parente | Sparse unmixing of hyperspectral data: The legacy of sunsal[END_REF] known as library mismatches. This problem was partially addressed in [START_REF] Iordache | Sparse unmixing of hyperspectral data[END_REF]. A bandwise scaling SY was suggested to apply to the dataset where the diagonal matrix S was obtained by solving

( X, Ŝ) = arg min X,S ||SY -DX|| F s.t. X ≥ 0, 1 T m X = 1 T n , (77)
A cyclic descent solution was suggested with the initialization of S = I turning the problem into iterative problems of least squares (w.r.t. S) and constrained least squares (w.r.t. X). The nonsparse fractional abundances of X and the corresponding spectral vectors of Y were removed after the first iteration to speed up the algorithm. Problem ( 77) is noncovex, therefore, the solution will be suboptimal. Additionally, the bandwise scaling S was estimated for a suboptimal sparse X.

In [START_REF] Fu | Semiblind hyperspectral unmixing in the presence of spectral library mismatches[END_REF], the library mismatch was incorporated into the sparse unmixing problem given by ( X, D ) = arg min

X 1 2 ||Y -D X|| 2 F + λ m i=1 ||x (i) || q 2,q , s.t. X ≥ 0, ||d i -d i || 2 < δ, i = 1, 2, ..., m. (78)
We should note that a surrogate function was suggested to be minimized instead of [START_REF] Ulyanov | Deep image prior[END_REF] to reduce the complexity of the algorithm (see [START_REF] Fu | Semiblind hyperspectral unmixing in the presence of spectral library mismatches[END_REF]). Additionally, the MUSIC pruning method was adapted (so-called robust MUSIC) to consider the mismatches given by arg min

z (d i -z) T P ⊥ (d i -z) ||d i -z|| 2 2 s.t. ||z|| 2 < δ. (79)
Recently, sparse unmixing using archetypal analysis (SUnAA) was proposed in [START_REF] Rasti | Sunaa: Sparse unmixing using archetypal analysis[END_REF], which addresses the spectral library mismatch from a different viewpoint. Instead of scaling data, SUnAA assumes that endmembers can be a convex combination of the library endmembers. Therefore, an extra matrix defines the contributions of the endmember from the library. The proposed linear model in [START_REF] Rasti | Sunaa: Sparse unmixing using archetypal analysis[END_REF] utilizes the advantages of the low-rank model and the sparse and redundant model defined as

Y = DBA + N, s.t. B ≥ 0, 1 T m B = 1 T r , and A ≥ 0, 1 T r A = 1 T n , (80)
where B ∈ R m×r determines the contributions of the endmembers from the the library D and A ∈ R r×n is the (lowrank) abundance matrix. To simultaneously estimate B and A SUnAA solves a minimization problem given by ( B, Â) = arg min

B,A 1 2 ||Y -DBA|| 2 F s.t.B ≥ 0, 1 T m B = 1 T r , andA ≥ 0, 1 T r A = 1 T n , (81)
SUnAA assumes that the unknown endmembers are a convex combination of the library's endmembers which leads to enforcing non-negativity and sum to one constraint on B. Unlike conventional sparse unmixing, the problem (81) is nonconvex.

A cyclic descent algorithm based on an active set method [START_REF] Chen | Fast and robust archetypal analysis for representation learning[END_REF] was proposed to solve [START_REF] Hong | Endmember-guided unmixing network (egu-net): A general deep learning framework for self-supervised hyperspectral unmixing[END_REF], iteratively. We should note that SUnAA is a parameter-free technique. Numerous techniques were tailored utilizing a combination of the aforementioned ideas, i.e., spatial regularization (TV), spatial/spectral weighted norm, collaborative norm, low-rank inducing, superpixel (segmentation), structured library/ group norm, and localization, library pruning [START_REF] Su | Superpixel-based weighted collaborative sparse regression and reweighted low-rank representation for hyperspectral image unmixing[END_REF]- [START_REF] Zhang | Spectral reweighting and spectral similarity weighting for sparse hyperspectral unmixing[END_REF]. Here, we have drawn the main ideas behind the majority of sparse techniques.

5) Shallow/Deep

Neural Networks: There are a few NNbased sparse unmixing approaches proposed in the literature. The pioneer method is called sparse unmixing using an unsupervised convolutional neural network (SUnCNN) [START_REF] Rasti | SUnCNN: Sparse unmixing using unsupervised convolutional neural network[END_REF]. As mentioned, the 1 penalties cannot be applied to the abundance while holding the ASC. This was addressed in [START_REF] Rasti | SUnCNN: Sparse unmixing using unsupervised convolutional neural network[END_REF] by implicitly applying an image prior and holding the ASC and ANC. It was shown that the problem of selecting a suitable prior for a sparse regression given by

X = arg min X 1 2 ||Y -DX|| 2 F + λφ(X) s.t.X ≥ 0, 1 T m X = 1 T n , (82)
could be moved to optimizing the parameters of a deep network (f θ ) given by

θ = arg min θ 1 2 ||Y -Df θ (Z)|| 2 F s.t. X = f θ (Z), (83)
where f θ is the deep CNN with parameters θ. and ANC and ASC can be enforced using a softmax layer. In [START_REF] Lin | Addressing spectral variability in hyperspectral unmixing with unsupervised neural networks[END_REF], an asymmetric encoder-decoder architecture is used for sparse unmixing. Instead of softmax, a sparse variation of softmax is used to avoid the full support of softmax while enforcing ASC. Spectral Variability Augmented Two-Stream Network (SVANT) was proposed in [START_REF] Zhang | Spectral variability augmented two-stream network for hyperspectral sparse unmixing[END_REF]. SVATN uses the convolutional encoder-decoder used in SUnCNN for both streams. The first stream uses a spectral library to estimate the abundances, and the second stream uses a spectral variability library to estimate spectral variation. A correlation-based variability extraction method was proposed to create a spectral variability library. A spectral-spatial feature learning encoder-decoder network was proposed in [START_REF] Kong | A global spectral-spatial feature learning network for semisupervised hyperspectral unmixing[END_REF]. The encoder utilizes two branches with spatial and spectral filters to extract the spatial and spatial features, respectively. The concatenated features go through three spectral-spatial residual blocks, which utilize MLP for spectral attention. On the other hand, the decoder is shallow and uses the endmember library to reconstruct data.

Algorithm unrolling-based approaches were also proposed for sparse unmixing. An unrolling-based shallow network was recently proposed for sparse unmixing in [START_REF] Qian | Deep unfolded iterative shrinkagethresholding model for hyperspectral unmixing[END_REF], [START_REF] Shao | Iviu-net: Implicit variable iterative unrolling network for hyperspectral sparse unmixing[END_REF]. An ADMM-based solution to the nonnegative 1 sparse regression problem (so-called SUnSAL), i.e., problem [START_REF] Bajjouk | Quantification of subpixel cover fractions using principal component analysis and a linear programming method: Application to the coastal zone of roscoff (france)[END_REF], was unrolled. An intermediate convolution was applied to the abundance to incorporate spatial information. A combination of SAD, MSE, and SID was used as the loss function to train the shallow network. In [START_REF] Qian | Spectral mixture model inspired network architectures for hyperspectral unmixing[END_REF], the iterative soft-thresholding algorithm (ISTA) [START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF] was unrolled to solve the nonnegative 1 sparse regression problem. SUnSAL was also unrolled in [START_REF] Zhou | Admm-based hyperspectral unmixing networks for abundance and endmember estimation[END_REF] for sparse unmixing.

VII. NONLINEAR UNMIXING APPROACHES

The linear mixture model is a simplified model and often fails in the case of intimate mixtures or/and when the light undergoes multiple reflections before reaching the sensor [START_REF] Bioucas-Dias | Hyperspectral remote sensing data analysis and future challenges[END_REF]. Alternatively, nonlinear models are used [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], [START_REF] Dobigeon | Nonlinear unmixing of hyperspectral images: Models and algorithms[END_REF].

Bilinear approaches are often used for double scattering. They assume that the received light at the sensor interacts with two materials. This is modeled using an extra mixing term by the Hadamard product between the endmembers of the materials. The Fan model [START_REF] Fan | Comparative study between a new nonlinear model and common linear model for analysing laboratory simulated-forest hyperspectral data[END_REF], [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF] is a variation of this model. One disadvantage of the Fan model is the generalization; therefore, it performs poorly for linearly mixed datasets. Polynomial post nonlinear mixing model (PPNM) [START_REF] Altmann | Supervised nonlinear spectral unmixing using a postnonlinear mixing model for hyperspectral imagery[END_REF], generalized bilinear model (GBM) [START_REF] Halimi | Nonlinear unmixing of hyperspectral images using a generalized bilinear model[END_REF], and linearquadratic model (LQM) [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF] were proposed to generalize the Fan model for the linear mixtures. However, they contain hyperparameters to describe the trade-off between the linear and nonlinear terms. In [START_REF] Gao | Using low-rank representation of abundance maps and nonnegative tensor factorization for hyperspectral nonlinear unmixing[END_REF], a nonlinear low-rank tensor unmixing algorithm was proposed to solve the GBM. Bilinear models have physical interpretation in some specific applications, such as canopy scenarios; however, they have several disadvantages. They have many parameters, and the estimated abundances are hard to interpret. Additionally, they often do not include self-interactions or consider the reflections from objects outside the instance field of view [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF]. Moreover, they are limited to secondary interactions. Therefore, several nonlinear mixing models have been developed, such as the multilinear mixing model (MLM) [START_REF] Heylen | A multilinear mixing model for nonlinear spectral unmixing[END_REF] and the p-linear (p > 2) mixture model (pLMM) [START_REF] Marinoni | A novel approach for efficient p-linear hyperspectral unmixing[END_REF]- [START_REF] Marinoni | Harmonic mixture modeling for efficient nonlinear hyperspectral unmixing[END_REF] were proposed for multiple interactions of the incident light.

Kernel methods have also been employed for nonlinear unmixing [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF]. They represent data in a higher dimensional space to linearize the problem. Hence, linear models become effective in high-dimensional space. In [START_REF] Broadwater | Kernel fully constrained least squares abundance estimates[END_REF], a kernel nonnegative matrix factorization (kernel-NMF) was proposed in [START_REF] Wu | A kernel spatial complexity-based nonlinear unmixing method of hyperspectral imagery[END_REF] for nonlinear unmixing. The FCLSU problem can be solved using support vector machines (SVMs). Pure pixels in data are the common support vectors that span the data simplex and allow the estimation of the abundances of all pixels enclosed within the simplex. Hence, Kernel SVMs were proposed for nonlinear spectral unmixing in [START_REF] Brown | Linear spectral mixture models and support vector machines for remote sensing[END_REF].

Radiative transfer models are mathematical physic-based mixing models that can be used to reconstruct intimately mixed materials' reflectance spectrum [START_REF] Heylen | A review of nonlinear hyperspectral unmixing methods[END_REF]. The inverse problems based on such models are very hard to solve. The Hapke model estimates the areal fractions of the materials within the mixture by transforming the reflectance spectra to their singlescattering albedos (SSA) and applying linear unmixing. The simplified version of the Hapke model was used for predicting intimate mixtures' composition [START_REF] Hapke | Bidirectional reflectance spectroscopy: 1. theory[END_REF], [START_REF] Hapke | The opposition effect of the moon: Coherent backscatter and shadow hiding[END_REF]. The proposed nonlinear unmixing methods based on Radiative transfer models often assume that spectral reflectance of the pure materials is available for estimating fractional abundances.

Deep networks were also proposed for nonlinear unmixing. Deep AE-based architectures have been broadly used for linear unmixing. In [START_REF] Wang | Nonlinear unmixing of hyperspectral data via deep autoencoder networks[END_REF], an AE network was proposed where the encoder utilizes an extra nonlinear layer to capture the nonlinearity of the data. A deep AE with multitask learning was suggested in [START_REF] Su | Deep autoencoders with multitask learning for bilinear hyperspectral unmixing[END_REF]. A long short-term memory-based autoencoder was proposed for PPNM in [START_REF] Zhao | Lstm-dnn based autoencoder network for nonlinear hyperspectral image unmixing[END_REF]. The proposed method in [START_REF] Zhao | Hyperspectral unmixing for additive nonlinear models with a 3-d-cnn autoencoder network[END_REF] exploits a 3-D CAE-based network for PPNM. A supervised AE was used in [START_REF] Shahid | Spatial-aware hyperspectral nonlinear unmixing autoencoder with endmember number estimation[END_REF] for Fan, bilinear, and PPNM. The Radial basis function (RBF) kernels and K-means clustering were used to estimate the number of endmembers and the endmember spectra, respectively. DL-based nonlinear unmixing techniques are often AE-based networks using PPNM, and as we discussed, they may have disadvantages of the bilinear models mentioned above [START_REF] Li | Model-based deep autoencoder networks for nonlinear hyperspectral unmixing[END_REF]- [START_REF] Lyu | Identifiability-guaranteed simplex-structured postnonlinear mixture learning via autoencoder[END_REF]. GAN architectures were also explored for nonlinear unmixing. In [START_REF] Tang | Hyperspectral nonlinear unmixing via generative adversarial network[END_REF], a cycle-consistent loss was used to ensure the reconstruction in addition to two GAN losses. A CNN was designed based on Hapke model (HapkeCNN) in [START_REF] Rasti | Hapkecnn: Blind nonlinear unmixing for intimate mixtures using hapke model and convolutional neural network[END_REF] to incorporate the physical model in the learning process.

VIII. EXPERIMENTAL RESULTS

We use three simulated and one real datasets. The simulated datasets were designed according to different mixing scenarios briefly explained in Table III. We avoid using the widely used benchmark datasets such as Samson and Jasper since their abundances are generated synthetically. For real-world experiments, we use the Cuprite dataset, a well-studied site with geological reference maps. The simulated experiments were carried out for 30 dB SNR. The tuning parameters were fine-tuned for the methods up to some levels. We should note that some methods have several parameters to be tuned; therefore, searching for the optimum is cumbersome. All the results are averaged over 10 experiments, and the standard deviations are shown by error bars. In experiments, we compare 20 unmixing methods from different categories as follows. For supervised methods, we used three endmember extraction/estimation techniques, i.e., VCA, SiVM, and SISAL with FCLSU and UnDIP. All six combinations of them were considered. For blind unmixing, we use PGMSU, MSNET, CNNAEU, ADMMNet, BayesianSMA, NMFQMV, MiSiC-Net, and EDAA. For sparse unmixing, we used SUnSAL, CLSUnSAL, MUA SLIC, S 2 WSU, SUnCNN, and SUnAA. The codes for all those methods were provided in HySUPP for reproducibility.

For the quantitative evaluation, we use the SRE in dB for estimated abundances given by SRE(X, X) = 20 log 10 X F X -X F .

A. Data Description 1) Synthetic Datasets with Spatial Structure: We simulated three data cubes (DC1, DC2, and DC3). DC1 was simulated using a linear mixing model with 5 endmembers selected from the USGS library and 75×75 pixels. The abundance maps are composed of five rows of square regions uniformly distributed over the spatial dimension. This dataset contains pure pixels for all endmembers. DC2 has 100×100 pixels and was simulated using a linear mixing model with 9 endmembers. The abundance maps were sampled from a Dirichlet distribution centered at a Gaussian random field to have piecewise smooth maps with steep transitions. Therefore, DC2 contains spectral variations. For DC1 and DC2, an endmember library D ∈ R 188×240 , composed of 240 spectral signatures were selected from the USGS library with a minimum pairspectra angle of 4.44°. DC3 was simulated with no pure pixels, and it has two mixed pixels on the facet of the data simplex. DC3 is composed of 105×105 pixels using the linear combination of six endmembers. For DC3, we use a library D ∈ R 188×498 composed of 498 spectral pixels from the USGS library. Note that we remove the water absorption and noisy bands, and therefore, the final pixels are of dimension p = 188.

2) Cuprite Dataset: The Cuprite dataset used in this paper contains 250× 191 pixels. Cuprite is a well-studied mineral site, and the dominant minerals are demonstrated using a geological ground reference therefore, the estimated abundance maps by different techniques can be compared visually. We use the same library as explained for DC3.

B. Experimental Results: Synthetic Datasets

Figure 13 demonstrates the unmixing results in terms of abundance SRE in DB for different techniques applied to DC1, DC2, and DC3. The outcomes of the results can be summarized as follows:

• For the pure pixel scenario, supervised techniques perform very well (Fig. 13,DC1). Overall, they perform better than the semisupervised and blind techniques. This confirms the importance of geometrical information for endmember extraction/ estimation techniques when there are pure pixels in the dataset. This is further confirmed in blind methods where MiSiCNet and NMFQMV, which exploit geometrical information, outperform the other blind techniques and provide competitive results compared to supervised ones. Sparse techniques show moderate results except SUnAA, which outperforms all the other techniques. The results confirm that the sparse unmixing techniques are not suitable when there are pure pixels for the endmembers in the dataset. We should mention that SUnAA does not match the characteristics of conventional sparse techniques. Even though SUnAA relies on a library, it uses a nonconvex optimization to endmembers. • On the other hand, for DC2, which contains spectral variations, as can be seen from Fig. 13, sparse unmixing techniques outperform the supervised and blind techniques. The results confirm that sparse unmixing techniques are more suitable for capturing the spectral variability. SUnAA outperforms the other technique. Note that SUnAA is a parameter-free technique. Supervised techniques outperform blind techniques due to the presence of pure pixels. • In the case of no pure pixel scenario (Fig. 13,DC3), blind unmixing techniques that exploit geometrical information outperform the other techniques. Sparse unmixing provides very poor results. Although SUnAA considerably significantly outperforms the performance of sparse unmixing techniques, it is very far from the best performance which is obtained by MiSiCNet. Among supervised techniques, SISAL (+ FCLSU/UnDIP) shows competitive results because it uses geometrical information to estimate the endmembers.

C. Experimental Results: Real Data

We selected three methods per category to conduct unmixing on the Cuprite dataset. Blind unmixing: MiS-iCNet, MSNet and NMFQMV. Semi-supervised: SUnAA, MUA SLIC, S 2 WSU. Supervised: UnDIP combined with SISAL, SiVM and VCA. We describe the hyperparameters that were fine-tuned for the following techniques. The hyperparameters set as default are omitted.

• MiSiCNet: λ = 100, projection=True.

• MSNet: α = 0.1, β = 0.1.

• MUA SLIC: β = 30, λ 1 = 0.001, λ 2 = 0.001, slic size = 200. • S 2 WSU: λ = 0.001.

• SISAL: τ = 1e-6. Fig. 14 (b) visually compares the estimated abundances for three dominant minerals, i.e., Chalcedony, Alunite, and Kaolinite. The comparison with the geological reference map Fig. 14 (a) reveals that the estimated abundances obtained by semisupervised methods show more resemblance to the reference map for all three minerals. SUnAA visually outperforms the other techniques, particularly in the case of Chalcedony. The blind unmixing methods can better estimate Chalcedony compared to MUA SLIC and S 2 WSU. This could be attributed to the mismatch of the endmember with the library's endmembers for this mineral. It is worth mentioning that SUnAA does not entirely rely on the library, and it estimates the endmember. Therefore, it can compensate for such a mismatch. More importantly, SUnAA is a parameter-free technique. We should note that selecting optimum parameters for the unmixing techniques is not a trivial task in real-world applications since the abundance RMSE cannot be computed.

IX. CONCLUSION, DISCUSSION, SUMMARY

Assuming that we have captured a spectral dataset and now have an unmixing problem in hand, we need to estimate the abundances of materials. The main question is which method to choose and which group of methods to select to tackle the problem. Indeed, the first step is to evaluate our problem and see if the linear mixing model or its variations are suitable for our problem. This decision needs prior knowledge of the physics of the problem. For instance, if you are dealing with intimate mixtures or close-range and microscopic scenarios, you should use nonlinear models. If you are dealing with macroscopic Earth observation problems, then linear models or their variations will be suitable. In some research, nonlinear models perform better than linear ones, however, one may pay attention to the selection of the model in a real-world application. Usually, linear models are more general.

Here, we clarify the keys to the success of each category. The success of supervised (or sequential) unmixing lies in the confidence of the endmember measured, extracted, or estimated (pure/no pure scenario). Therefore, we should not use the supervised method if we are not confident about the endmembers. In other words, supervised methods perhaps are the best choice for endmembers with high confidence. When we have prior information on the material in the scene and a welldesigned endmember library, semisupervised unmixing could be successful. Semisupervised are also suitable to capture the spectral variability. The success of semisupervised unmixing lies in the quality of the endmember library. Blind unmixing methods should be selected when there is no library, no pure pixels in the data set (including highly mixed scenarios), or the confidence of the measure, selected, or estimated endmembers Abundances SRE (dB) Abundances SRE (dB) Abundances SRE (dB) is low. They should be used with caution, and the estimated endmembers should always go through physical interpretation.

Despite the considerable advances, spectral unmixing can be considered one of the most challenging tasks in hyperspectral analysis. Here, we describe some of the main unmixing challenges.

• The linear models are much more general than the nonlinear models. However, from one dataset to another dataset, their performance may significantly decrease. • Parameter selection considerably affects the performance of the unmixing methods. The selection of optimum parameters is very challenging when it comes to realworld datasets. • The performance of the linear unmixing methods often degrades by increasing the number of endmembers. The linear unmixing may fail for a dataset with a higher number of endmembers. • Spectral variability is still a big challenge and may considerably downgrade the performance of linear unmixing. • One of the main challenges is the absence of a real dataset with ground truth. • Multitemporal and multisource spectral unmixing are also challenging tasks. • Considering the large volume of hyperspectral data, scalable unmixing is the main key for global monitoring.

Fig. 1 .

 1 Fig. 1. Comparisons of mixed and pure pixels in hyperspectral data.

Fig. 2 .Fig. 3 .

 23 Fig. 2. Macroscopic versus microscopic assumption. This figure illustrates three major assumptions in hyperspectral imagery. (a) The linear assumption is when the light interacts with the materials only once. (b) The bilinear assumption when the light interacts with a maximum of two materials (c) intimate mixture when the light interacts with more than two materials.

Fig. 4 .Fig. 5 .

 45 Fig. 4. Publications over time based on IEEE Xplore keyword search tool using "Hyperspectral Unmixing" as input.

Fig. 6 .

 6 Fig. 6. Schematic of sensing and linear unmixing. (a) Sensing a mixed pixel (b) Reflectance of the mixed pixel and (c) Schematic of the linear unmixing.

Fig. 7 .

 7 Fig. 7. Graphical representation of different types of linear unmixing. Supervised and blind unmixing use the same low-rank mixing model, while sparse unmixing uses the sparse and redundant linear model. Blind unmixing estimates both endmembers and abundances simultaneously.

Fig. 9 .

 9 Fig. 9. Taxonomy of the linear unmixing methods described in this paper. The top categorization is based on the prior knowledge of endmembers.

Fig. 11 .

 11 Fig. 11. Pure pixel versus no pure pixel assumption. (a) Pure pixels scenario (b) No pure pixels scenario with enough pixels on the facets of the simplex (c) endmember bundle, and (d) No pure pixels scenario without pixels on the facets of the simplex (highly mixed scenario).

Fig. 12 .

 12 Fig. 12. The intuition of how the three geometrical regularizers enforce the endmembers using (a) TV, (b) center, and (c) boundary terms to change an initial simplex (blue) towards the endmember simplex (red). The black arrows show the direction of the force.

C

  N N A E U E D A A M S N e t P G M S U A D M M N e t B a y e s ia n S M A N M F Q M V M iS iC N e t S U n S A L C L S U n S A L M U A _ S L IC S 2 W S U S U n C N N S U n A A S IS A L + F C L S S IS A L + U n D IP S iV M + F C L S V C A + F C L S V C A + U n D IP S iV M + U n

  M N e t E D A A B a y e s ia n S MA N M F Q M V M S N e t C L S U n S A L S U n S A L M U A _ S L IC S U n C N N S 2 W S U S U n A A S IS A L + F C L S S IS A L + U n D IP S iV M + F C L S S iV M + U n D IP V C A + F C L S V C A + U n

  Abundances SRE on DC2 (SNR=30dB) blind semi-supervised supervisedC N N A E U A D M M N e t P G M S U B a y e s ia n S M A M S N e t E D A A N M F Q M V M iS iC N e t S U n S A L C L S U n S A L S 2 W S U M U A _ S L IC S U n C N N S U n A A V C A + U n D IP V C A + F C L S S iV M + U n D IP S iV M + F C L S S IS A L + U n D IP S IS A L + F C L S

Fig. 13 . 105 )

 13105 Fig.[START_REF] Zare | Endmember variability in hyperspectral analysis: Addressing spectral variability during spectral unmixing[END_REF]. Comparing abundance SRE (↑) in dB using different unmixing techniques applied to (from top to bottom) synthetic DC1, DC2, and DC3.

Fig. 14 .

 14 Fig. 14. Estimate abundances obtained by applying different unmixing techniques to Cuprite compared with the geological reference map.

TABLE I THE

 I LIST OF UNMIXING METHODS INCLUDED IN HYSUPP WITH THEIR CORRESPONDING DEPENDENCIES. WE HIGHLIGHTED OUR MAIN CONTRIBUTIONS AND THE LINK TO THE ORIGINAL IMPLEMENTATIONS.

	Method	Original implementation	Python	Dependencies	GPU	Contributions
	FCLSU [21]	pysptools		numpy, cvxopt

TABLE II SELECTION

 II OF G AND O FOR THE UNIFIED GEOMETRICAL REGULARIZER.

ACKNOWLEDGMENTS

AZ, JM, and JC were supported by ANR 3IA MIAI@Grenoble Alpes (ANR-19-P3IA-0003). AZ and JM were supported by the ERC grant number 101087696 (APHELAIA project).

APPENDIX MINIMUM VOLUME SIMPLEX USING TOTAL VARIATION

The total variation penalty enforces the data simplex to have minimum volume by pulling the endmembers towards each other. Assuming ē = 1 r e1 r we have

where the third term is zero. Therefore,