Archetypal Analysis for Sparse Representation-based Hyperspectral Sub-pixel Quantification

02/08/2018 ∙ by Lukas Drees, et al. ∙ University of Bonn 0

The estimation of land cover fractions from remote sensing images is a frequently used indicator of the environmental quality. This paper focuses on the quantification of land cover fractions in an urban area of Berlin, Germany, using simulated hyperspectral EnMAP data with a spatial resolution of 30m×30m. We use constrained sparse representation, where each pixel with unknown surface characteristics is expressed by a weighted linear combination of elementary spectra with known land cover class. We automatically determine the elementary spectra from image reference data using archetypal analysis by simplex volume maximization, and combine it with reversible jump Markov chain Monte Carlo method. In our experiments, the estimation of the automatically derived elementary spectra is compared to the estimation obtained by a manually designed spectral library by means of reconstruction error, mean absolute error of the fraction estimates, sum of fractions, R^2, and the number of used elementary spectra. The experiments show that a collection of archetypes can be an adequate and efficient alternative to the manually designed spectral library with respect to the mentioned criteria.



There are no comments yet.


page 3

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

I Introduction

The estimation of the degree of imperviousness as an indicator of environmental quality is subject of current research towards a time and cost efficient monitoring of urban areas [1]. Due to increasing land consumption in cities in the recent years, which has negative effects on the natural water cycle, the monitoring of land use in those areas is important [2].

Remote sensing data, such as imaging spectroscopy, builds a valuable basis to comprehensively map urban areas and quantify the imperviousness based on the spectral information (e.g., [3, 4]). Especially, hyperspectral imagery is a suitable source for mapping of such areas, because it offers a high spectral separability of different materials. However, generally, the temporal and spatial resolution is limited in comparison to sensors with lower spectral resolution. These limitations are partially overcome with the launch of missions such as Environmental Mapping and Analysis Program (EnMAP), which increases the availability of hyperspectral data and the temporal resolution [5]. Nevertheless, due to its spatial resolution, the provided data is mainly characterized by spectrally mixed pixels, which demands sophisticated sub-pixel quantification approaches in order to estimate the fraction of various land cover classes in each pixel.

In this context, several approaches have been developed comprising regression approaches [6, 7], probabilistic classification methods [8, 9, 10], and the usage of spectral libraries for spectral mixture analysis [11, 12]. An overview of a wide variety of unmixing approaches can, for example, be found in [13]. While the latter approach needs a spectral library containing elementary spectra of known materials, the first two approaches also require mixed spectra for learning an appropriate model. These mixed spectra can be derived from the image using information about known mixed pixels, or from synthetically mixed pixel, as it has been presented in [8] and [6].

When using spectral libraries, a critical step is the extraction of the elementary spectra. A manual extraction is time-consuming and requires human expert-knowledge and therefore, automatic extraction techniques have been an active field of research during the past decade (e.g., [13, 14]). Most of the algorithms rely on the assumption that the elementary spectra lie on a convex hull or a convex polytope enclosing the data distribution (e.g., [15, 16]). Based on this assumption all data samples can be reconstructed by a non-negative linear combination of the elementary spectra. A promising approach from this group is the so-called archetypal analysis, which searches for extreme points (also known as archetypes) in the data distribution (e.g., [17, 18, 19, 20]). Archetypal analysis has already been successfully applied in the field of sport analytics [21], plant phenotyping [22] or text analysis [23]. A valuable extension to archetypal analysis is presented by [24], in which extreme points are extracted in the kernel space, enabling an efficient nonlinear unmixing. Besides the actual determination of elementary spectra, other challenges exist which need to be tackled: The number of elementary spectra is unknown beforehand and thus, a suitable amount of spectra needs to be extracted to make the set representative enough, but also small enough to keep the sub-pixel quantification robust and efficient. Moreover, many extraction techniques depend on the initialization and thus, a strategy needs to be defined to ensure a stable result (e.g., [25, 26]).


Part of HyMap image data visualized as RGB-image with the wavelengths R =

nm, G = nm and B = nm
(b) Reference image with four classes impervious surface, vegetation, soil and water
(c) Simulated EnMAP data visualized with the wavelengths as specified above (bands: RGB = 11,5,1)
Fig. 1: Berlin-Urban-Gradient dataset

In this paper, we address the challenge of automatically finding a representative set of elementary spectra, including the automatic determination of the number of elements, for sub-pixel quantification. We perform the sub-pixel quantification using a freely available simulated EnMAP scene of an urban area in Berlin, Germany [27]111, aiming at the estimation of a fraction map containing the classes impervious surface, vegetation, soil and water. In order to determine the class fractions, we use sparse representation with non-negativity, L-sparsity and sum-to-one constraint. We exploit archetypal analysis to extract elementary spectra in a fully automatic and unsupervised way. Moreover, we perform archetypal analysis by simplex volume maximization (SiVM), which states an efficient selection method [28, 29].

The main contribution of this work is the combination of archetypal analysis with reversible jump Markov chain Monte Carlo method (rjMCMC, [30]) to obtain a spectral library of high representational power, yet a small number of elementary spectra. Using this approach, we are able to determine the best set of spectral library elements regarding a pre-defined criteria, as well as the number of elements. Moreover, the approach of [24] is applied to select the archetypes in kernel space, resulting in archetypes lying on the concave hull of the data distribution in the original feature space. Our experiments confirm that these archetypes are more suitable for sub-pixel quantification than archetypes extracted from the convex hull exclusively. To illustrate the usefulness of our proposed approach, we analyze various kinds of automatically derived spectral library and compare them to a manually designed spectral library. Our presented approach is flexible regarding the chosen estimation technique, such that the constrained sparse representation can be replaced by other approaches commonly used in the unmixing community ([13]), or other constraints such as sparsity induced by L-norm ([31]). As such, archetypal analysis can be replaced by, e.g., endmember extraction methods presented in [13], and deliver the input for the rjMCMC method.

Ii Data

Our studies are performed using the Berlin-Urban-Gradient dataset [27], illustrated in Fig. 1. The dataset consists of two hyperspectral images of different spatial resolution, two simulated EnMAP scenes of different spectral resolution, a manually designed spectral library, reference land cover information and reference fractions for evaluation, which are explained in more detail in the following paragraphs.

Fig. 2: Significant spectra of the spectral library (Lib). Asphalt and red clay tile roof in the class imp. surface, grass and tree for vegetation, bare soil and natural water.
Imp. surface Vegetation Soil Water
39 31 4 1

TABLE I: Composition of 75 spectra in the manually designed spectral library (LibCom)

Ii-1 HyMap Data

The dataset contains two hyperspectral images, one with a spatial resolution of m and one with a spatial resolution of m, whereby we use the higher resolution image in our work. Both images were acquired with the Hyperspectral Mapper (HyMap), and serve as basis for the extraction of the archetypes. It offers a high variety of urban land use and land cover patterns in the study site Southwest of Berlin, Germany. The external conditions were a cloudless sky at solar noon and a minimal possible altitude (maximal altitude according to the second lower resolution image). The HyMap image consists of 126 spectral bands, however, 15 noisy bands were removed resulting in 111 bands used for this work. Moreover, other preprocessing steps were performed, as encompassed system correction [32], atmospheric correction and parametric geocoding [33]. The observed wavelengths range from 0.45m to 2.5m, showing a high spectral information diversity, which enables a detailed analysis of the urban structure. In addition to the m-resolution HyMap image, a manually derived reference image containing reference pixel is provided with valuable spectra, each labeled with one of the four land cover classes: impervious surface, vegetation, soil and water. The reference information is manually obtained by using digital orthophotos and cadastral data. For our experiments, the set of elementary spectra is extracted from this data to obtain the class membership of the elementary spectra, however, the extraction can be performed on unlabeled data with limited human user interaction.

Ii-2 Simulated EnMAP Data

EnMAP is German hyperspectral satellite mission, which will start not earlier than 2018, with a focus on Earth environmental observations in a global scale. Based on the HyMap data, an EnMAP scene was simulated using EnMAP end-to-end simulation tool (EeteS, [34]) with two different spectral resolutions (111 and 244 bands). Just as the HyMap data, both EnMAP images have a spectral resolution ranging from m to m, where we use bands for a better comparability. The spatial resolution of m is lower than the resolution of the HyMap scene, and thus the mixing of land cover classes are more apparent. For evaluation purposes, EnMAP pixels were obtained from the simulation tool, containing the fractions of the land cover classes ranging from to .

Ii-3 Manually Designed Spectral Library

The spectral library is a manually designed collection of pure spectra obtained from the HyMap image. The spectra contain different land cover classes with impervious surface spectra (different types of roof, pavement, tartan, pool water), vegetation spectra (grass, tree), soil spectra (uncovered ground, sand) and one natural water spectrum. All together, it is a balanced library for urban structures with impervious and pervious spectra. All spectra in the library can be assigned to a hierarchical urban classification scheme, which was developed by [35]. First, an initial collection of spectra was selected by expert knowledge. Afterwards a two-step filtering was performed, in which the variability between the spectra is maximized in consideration of spectral variability of materials, brightness and shading effects. Moreover, in an iterative process those spectra for the final subset were selected that best describes the spectral diversity in a specific neighborhood. More details can be found in [27].

Iii Methods

The following section describes the sub-pixel quantification using sparse representation and archetypal analysis. Archetypal analysis determines the extreme points of the data distribution, which are used within the sparse representation approach to estimate the fractions of land cover classes in each pixel. We have given a -dimensional data matrix X, in which is the number of -dimensional reference pixels with given land cover class . Moreover, we have given a test set of data samples with known land cover fractions for evaluation purposes, as presented in Sec. II-2.

Iii-a Sparse Representation

In order to determine the sub-pixel fractions, we use sparse representation with non-negative least squares optimization. In terms of sparse coding, a sample is represented by a weighted linear combination of a few elements taken from a -dimensional dictionary D, such that with being the reconstruction error. The dictionary

contains elementary spectra, such that this approach is identical to the linear mixing model. The coefficient vector, comprising the activations, is given by

. The activations are interpreted as class fractions for sub-pixel quantification. The optimization problem for the determination of optimal is given by


where the terms in (2) are the non-negativity constraint, the sum-to-one constraint, and the sparsity constraint. Generally, non-negativity alone leads to a sparse solution, however, in order to ensure a strict fulfillment of the sparsity constraint we use a backward selection procedure in which the activations with the smallest values are set to zero in a greedy manner.

Iii-B Archetypal Analysis

Archetypal analysis is a suitable method to determine the elements of D, where each archetype serves as one dictionary element. The extraction of the archetypes, collected in a -dimensional matrix of , , is carried out by SiVM, which is an efficient method to determine the archetypes of the data distribution. This approach is aiming on an approximation of the convex hull, where all archetypes are located on it.

In order to find the first archetype, the approach is initialized with a random or pre-defined vector . The pixel with maximal distance to , is defined as first archetype , defined by


with being the Euclidean distance function between the spectral features of the archetype and the pixel . Further archetypes are specified sequentially, such that the volume of the simplex becomes maximized with each additional archetype. Since the volume operation is too computational intense, instead the archetypes are selected to have maximum distance to all previously detected ones, using


The stopping criterion is generally chosen to be the number of archetypes.

We further use the approach of [24] and transform X

into kernel space using a Gaussian radial basis function kernel with a hyperparameter

describing the width of the Gaussian kernel. In this way, archetypes are selected which lie on the the concave hull rather than the convex hull.

The disadvantage of archetypal analysis is that the final set depends on the initialization point, and as a result, there is no unique solution to the final set. Especially, if the number of archetypes in the dictionary is low, various solutions lead to significantly different accuracies. Moreover, the number of archetypes is generally not know beforehand, and depends on the number of informative dimensions and the variability of the data.

Iii-C Reversible Jump Markov Chain Monte Carlo Method

To overcome this problem, we propose to use an optimization procedure to find the best set of archetypes from a large set of pre-selected ones, called the initial set. Under the assumption that a suitable archetypal set is able to represent the test set with a low reconstruction error, our task is to find the set of archetypes which minimizes the energy


where is obtained by using the current set of archetypes as dictionary D. The energy U is a complex function with rough landscape and unknown dimensionality due to the unknown number of archetypes. Therefore, we optimize with rjMCMC coupled with simulated annealing to find the global optimum. Introducing the temperature parameter , the optimizer is given by


While MCMC is dedicated to sample from complex functions, simulated annealing allows to make a point estimate of its global optimum. Using simulated annealing we create a Markov chain, such that the samples explore the whole state space in the beginning and gradually concentrate around the global optimum of the energy function U. In this way we avoid trapping into local optima, as it is usually the case for greedy algorithms. We use the so-called birth an death algorithm [36] to sample from a restricted sample space of possible sets of archetypes, which turns out to be a special type of Green’s rjMCMC sampler [30]. The sample space is restricted by a Poisson prior on the expected number of archetypes, as presented in [37]. We further introduce an upper bound on the the number of selected archetypes.

Iv Experiments

Iv-a Experimental Setup

In our experiments, the simulated EnMAP data is reconstructed by sparse representation with the before mentioned constraints using the following libraries:

  1. the manually designed spectral library (LibCom),

  2. an optimized set of the manually designed spectral library LibRed (using rjMCMC),

  3. a library containing extracted archetypes in the original feature space initialized by the mean vector (AA-M-Lin),

  4. a library containing extracted archetypes in the kernel space initialized by the mean vector (AA-M),

  5. a library containing extracted archetypes in the kernel space initialized by a random vector (AA-Rand),

  6. a library containing an accumulated set of multiple single sets of archetypes, obtained in the kernel space with random initializations (AA-Full), and

  7. an optimized set of AA-Full denoted by AA-Opt (using rjMCMC).

We compare our results to the regression and classification methods presented in [8]

, namely support vector machine (

SVM), import vector machine (IVM), support vector regression (SVR), and multi-output support vector regression (MSVR

). The aim of the experiments is to show the suitability of the spectral libraries for sub-pixel quantification. Moreover, the goal is to show that the set of automatically derived spectral libraries using archetypal analysis and rjMCMC achieve similar results than the manually designed spectral library. The sub-pixel quantification is evaluated with the given reference mixing fractions. Additionally, the number of used archetypes is evaluated and discussed. We run all experiments including a random component 10 times and report the average result and standard deviation.

As pre-processing step, the dataset is outlier-cleaned using the local outlier factor approach

[38], using

neighbors and a quantile of

as threshold on the pairwise Euclidean distances. The width of the Gaussian kernel is chosen to be of the average standard deviation over all bands. In order to choose a suitable sparsity value , we use Elbow method to analyze the reconstruction error of the test set obtained by the manually designed library as well as the automatically derived library containing archetypes. Fig. 3 shows that the reconstruction error will not significantly decrease for , and thus we choose this value as upper bound for the sparsity constraint. For rjMCMC we choose the prior on the expected number of archetypes to be when using all labeled data, and when using the manually designed spectral library. The upper limit of archetypes is chosen to be .

Fig. 3: Reconstruction error evaluated on the test set for different sparsity values .

Iv-B Evaluation

We use several statistics for evaluation. First, the reconstruction error of the sparse representation is provided, indicating the ability of the elementary spectra to represent the EnMAP pixels. Moreover, the result of the sparse representation is evaluated by a comparison between reference and actual estimated fraction coefficients. This is done by means of the overall and class-wise mean absolute error () between reference and actual coefficients.


with being the reference fractions, being the estimated fractions, and being the number of evaluated pixels. We obtain by summing up all coefficients belonging to elementary spectra of the same class. The smaller the , the better the reconstruction. Moreover, the root mean square error () is provided, which is defined by


Besides this, we provide the coefficient of determination () for each class. It is calculated as the squared correlation coefficient between and .

Iv-C Results and Discussion

Data Approach Amount Class-wise MAE [%] Ø MAE [%] Ø RMSE [%]
of elements Imp. surface Vegetation Soil Water
SVM 75 14.81 (15.22) 16.97 (15.97) 04.84 (19.01) 04.25 (66.33) 10.22 17.36
IVM 75 15.21 (15.91) 15.31 (15.59) 02.09 (20.98) 01.98 (35.51) 08.65 15.63
SVR 75 11.73 (11.51) 12.20 (11.92) 08.36 (12.78) 07.95 (07.68) 10.06 13.74
MSVR 75 11.33 (11.45) 10.99 (11.49) 02.17 (13.45) 03.19 (09.70) 06.92 11.32
LibCom 75 13.00 (13.79) 09.61 (09.87) 02.01 (14.37) 07.96(07.88) 08.15 12.66
LibRed 44.4 14.00 (14.86) 10.58 (10.85) 02.23 (15.07) 07.46 (09.07) 08.57 13.68
9.5 1.45 (1.57) 0.35 (0.39) 0.64 (0.61) 2.00 (1.71) 0.96 1.22
X AA-M-Lin 40 25.09 (26.23) 15.05 (15.70) 01.62 (20.38) 33.00 (09.65) 18.69 23.65
X AA-M 75 25.50 (27.00) 22.44 (22.93) 01.80 (02.04) 08.87 (20.39) 14.65 20.01
X AA-Rand 75 21.05(21.58) 17.86 (18.27) 01.77 (15.69) 08.52 (26.29) 12.30 18.17
X AA-Full 142 49.06 (52.71) 32.50 (33.20) 02.38 (26.16) 2.84 (68.13) 21.70 31.45
X AA-Opt 81.7 16.29 (16.21) 14.21 (14.70) 03.14 (24.09) 07.00 (09.73) 10.16 16.08
13.6 0.69 (0.76) 0.65 (0.60) 0.60 (2.91) 1.50 (1.38) 0.51 0.47

TABLE II: Evaluation results of EnMAP sub-pixel quantification. Overall (Ø) and class-wise (in brackets computed without zero-reference fractions) and overall . Sparse representation-based estimations with various obtained spectral libraries (explanations see Sec. IV-A) are compared to regression (SVR, MSVR) and classification methods (IVM, SVM). The standard deviation values are indicated by .
Fig. 4: Scatter plots obtained from LibCom (top row), LibRed (middle row) and AA-Opt (bottom row) representing the class-wise fraction estimates of impervious surface, vegetation, soil and water, opposed to the reference fractions. A correct estimated pixel lies on the dashed line. The solid line represents the respective least-square regression line for the scatter points. The formula for the regression line the and coefficient of determination () is also given in each plot.
(a) Full manually designed spectral library LibCom
(b) Reduced manually designed spectral library LibRed
Fig. 5: The frequency of the usage of elementary spectra for the reconstruction of the reference EnMAP pixels. Different colors indicate the class membership of the reference pixels (I: Imp. surface, V: Vegetation, S: Soil, W: Water). The maximal possible height of a bar is , if the elementary spectrum is involved in the reconstruction of all reference EnMAP pixels. The lower bright part of each bar specify the sum of the fractions (i.e., the estimated coefficients) over all pixel.

Table II presents the results of the EnMAP sub-pixel quantification. In the first block, the results of [8] are shown for comparability with regression and classification results. The middle block shows the results obtained by the full and reduced manually designed spectral library. The bottom part of the table presents the results obtained by archetypal analysis.

Overall, the values show that all spectral libraries which are obtained from the full manually designed library achieve similar and satisfactory results. It can be seen that in comparison to the results of the regression and classification approaches presented in [8], sparse representation with the mentioned constraints show equivalent results. Both spectral libraries have a small average , which is also obtained by IVM and MSVR. Remarkably, LibCom and LibRed show similar results, however, with LibRed having much fewer elements. LibCom has slightly lower class-wise for the classes impervious surface, vegetation and soil, where LibRed obtains better results for water. This indicates that the presented rjMCMC method is able to find a suitable subset of elementary spectra having nearly the same representational power than the manually designed library with all elements. Moreover, it can be stated from this result that the full spectral library contains redundant elementary spectra which can be discarded for sub-pixel quantification.

These findings are underlined by Fig. 4, which shows the scatter plots representing the class-wise fraction estimates opposed to the reference fractions. It can be observed that soil and water have a high proportion of reference fractions and lesser > fractions than impervious surface and vegetation. The impervious surface scatter plot show that for LibCom the two lines intersect nearly at the 15% point on the x-axis. This means that the estimated fractions are usually too small for all pixels which have a degree of impervious surfaces over , and otherwise too high for smaller values than the intersection point. In comparison to LibRed, the reference values are also mostly underestimated, but slightly better than LibCom. The soil scatter plots have large discrepancies between the true 1:1 line (dashed line) and the estimated least-squares regression line (solid line). Besides this, there is a high amount of fractions, which are estimated with up to when using the LibCom spectra and fractions when using LibRed spectra. Furthermore, all pixels with soil-fractions are clearly underestimated equally. This is also the case for the regression and classification approaches presented in [8]. We assume that the small number of elementary spectra is insufficient for the spectral diversity of soil pixels. Finally, the described observation with high estimated fractions also occurs in the Water class with values over .

Figure 5 illustrates the frequency of the usage of each spectrum for reconstruction in order to determine the fractions of EnMAP pixels. The results obtained by LibCom are shown on the left and the results for LibRed on the right side. Moreover, the lower bright part of each bar indicates the sum over all fractions from which the total proportion of each land cover class in the study area can be derived. LibCom shows that there are significant elementary spectra, e.g., 18, 19 (impervious surface), 45, 48, 52, 54 (vegetation) and 75 (water), which show a high fraction in the reconstruction. In the other side, some elementary spectra in this plot are infrequently used. It can be observed that there is generally no dependence between the height of a bar, i.e. the number of non-zero estimated fractions for the spectra, and its overall fraction’s sum for all pixels. A striking example is spectra 19. Spectrum 19 in the LibCom library is used for more than pixel reconstructions, but its fraction’s sum is lower in comparison to, e.g., spectrum 18, which is used more infrequently. However, LibRed has a more balanced composition of elementary spectra with fewer bar heights, which are close to zero. Especially noticeable are the differences in the bar for the class water between LibCom and LibRed. A water spectrum has the special characteristic that there are only small reflectance values over all bands, so that it acts almost as a linear factor. Because of this, it is well suited to support all reconstructions, resulting in a high overall fractions’ sum. Its value in the LibCom is high with a share in over reconstructions, in comparison to LibRed which shows a bar height of just over 600 and a lower fractions’ sum. We assume the reason for this is the smaller number of elementary spectra, where oftentimes a reconstruction with water spectra is helpful regardless of the presence of water in the pixel.

The results in Tab. II obtained from various spectral libraries based on the labeled data set X are more variable, indicating their dependency of a proper choice of hyperparameters such as the number of archetypes and the selection approach. For example, the number of selected archetypes in set AA-M-Lin is , which is smaller in comparison to AA-M. Both sets are initialized with the mean vector of all samples in X, where the first set is extracted in the original feature space and the latter one in the kernel space. We observed that a selection of more than archetypes in the set AA-M-Lin means that some archetypes lie in the center of the dataset X, since they have the maximal distance to the previously selected archetypes. The set AA-M-Lin consists of impervious surface spectra, vegetation spectra, soil spectra and one water spectrum. The total number of selected archetypes in AA-M is chosen to be , as for the manually designed library LibCom. For all classes except the water class, the results obtained by AA-M are worse than these ones obtained by AA-M-Lin. Also the results for the randomly initialized archetypal sets AA-Rand show only slightly better results. Moreover, we observed that oftentimes the selected archetypal set with random initialization does not contain all land cover classes. Thus, as indicated by the results, a single set on selected archetypes is not suitable enough for an accurate sub-pixel quantification. Therefore, in order to create a higher diversity of archetypes, various initializations are used for SiVM and accumulated into a larger set. The set AA-Full contains different archetypes, however, the high amount of elementary spectra results in high values. The best result based on X is obtained by AA-Opt, which is the reduced AA-Full set. This set has the best overall and the class-wise are in most cases comparable to the approaches based on . The final number of elements in set AA-Opt is on average, which is similar to the number of elements in the manually designed spectral library.

Fig. 4 underlines these findings. For impervious surface, vegetation and water the results are similar for all three libraries. However, also for soil the set AA-Opt has large discrepancy between the true 1:1 line (dashed line) and the estimated least-squares regression line (solid line), which is more distinct than for LibCom and LibRed.

V Conclusion

The quantification of sub-pixel fractions in remote sensing images is a relevant task to assess the environmental quality, for example, in urban areas. This paper presents a sub-pixel quantification of the urban area of Berlin, Germany, into four land cover classes impervious surface, vegetation, soil and water, using a manually designed spectral library and various kinds of automatically derived spectral libraries. We perform the estimation of the fractions by using sparse representation with non-negativity, sparsity, and sum-to-one constraints. We use archetypal analysis by simplex volume maximization to automatically derive the elements for a library, and apply reversible jump Markov chain Monte Carlo method to find a small, yet representative set of suitable elementary spectra. The archetypes are extracted from HyMap data with given reference information for all land cover classes. As our experiments suggest, the extracted archetypes are suitable to serve as spectral library for sub-pixel quantification. Moreover, in contrast to a manually designed library, the automatically derived spectral library can be easily extracted with no or limited human user interaction, and the library is specifically adapted to the current image characteristics. In case no reference data is given, the interpretation can be done in a fast way by assigning land cover classes to the extracted archetypes using expert knowledge. The presented approach is flexible regarding the chosen estimation technique for the sub-pixel fractions, and can also be combined with spatial information, which may further increase the approximation ability and accuracy of the fraction estimates.


The authors would like to thank the reviewers for their valuable comments, and Akpona Okujeni, Sebastian van der Linden and Patrick Hostert for providing the dataset. Moreover, the authors would like to thank Johannes Rosentreter and Andres Milioto for valuable discussions.


  • [1] Q. Weng, “Remote sensing of impervious surfaces in the urban areas: Requirements, methods, and trends,” Remote Sensing of Environment, vol. 117, pp. 34–49, 2012.
  • [2] G. Wessolek, “Bodenüberformung und -versiegelung,” Handbuch der Bodenkunde, 2001.
  • [3] X. Huang, Q. Lu, and L. Zhang, “A multi-index learning approach for classification of high-resolution remotely sensed images over urban areas,” ISPRS Journal of Photogrammetry and Remote Sensing, vol. 90, pp. 36–48, 2014.
  • [4] S. Roessner, K. Segl, M. Bochow, U. Heiden, W. Heldens, and H. Kaufmann, “Potential of hyperspectral remote sensing for analyzing the urban environment,” Urban remote sensing: Monitoring, synthesis and modeling in the urban environment, pp. 49–61, 2011.
  • [5] L. Guanter, H. Kaufmann, K. Segl, S. Foerster, C. Rogass, S. Chabrillat, T. Kuester, A. Hollstein, G. Rossner, C. Chlebek et al., “The enmap spaceborne imaging spectroscopy mission for earth observation,” Remote Sensing, vol. 7, no. 7, pp. 8830–8857, 2015.
  • [6] A. Okujeni, S. van der Linden, S. Suess, and P. Hostert, “Ensemble learning from synthetically mixed training data for quantifying urban land cover with support vector regression,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2016.
  • [7] F. Priem, A. Okujeni, S. van der Linden, and F. Canters, “Use of multispectral satellite imagery and hyperspectral endmember libraries for urban land cover mapping at the metropolitan scale,” in SPIE Remote Sensing.   International Society for Optics and Photonics, 2016.
  • [8] J. Rosentreter, R. Hagensieker, A. Okujeni, R. Roscher, and B. Waske, “Sub-pixel mapping of urban areas using enmap data and multioutput support vector regression,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2017, to appear.
  • [9] S. Suess, S. van der Linden, P. J. Leitao, A. Okujeni, B. Waske, and P. Hostert, “Import vector machines for quantitative analysis of hyperspectral data,” IEEE Geoscience and Remote Sensing Letters, vol. 11, no. 2, pp. 449–453, 2014.
  • [10]

    J. Zhu and T. Hastie, “Kernel logistic regression and the import vector machine,”

    Journal of Computational and Graphical Statistics, 2012.
  • [11] B. Somers, G. P. Asner, L. Tits, and P. Coppin, “Endmember variability in spectral mixture analysis: A review,” Remote Sensing of Environment, vol. 115, no. 7, pp. 1603–1616, 2011.
  • [12] R. L. Powell, D. A. Roberts, P. E. Dennison, and L. L. Hess, “Sub-pixel mapping of urban land cover using multiple endmember spectral mixture analysis: Manaus, brazil,” Remote Sensing of Environment, vol. 106, no. 2, pp. 253–267, 2007.
  • [13] J. M. Bioucas-Dias, A. Plaza, N. Dobigeon, M. Parente, Q. Du, P. Gader, and J. Chanussot, “Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 5, no. 2, pp. 354–379, 2012.
  • [14] M. A. Veganzones and M. Grana, “Endmember extraction methods: A short review,” in International Conference on Knowledge-Based and Intelligent Information and Engineering Systems.   Springer, 2008, pp. 400–407.
  • [15] T.-H. Chan, W.-K. Ma, A. Ambikapathi, and C.-Y. Chi, “A simplex volume maximization framework for hyperspectral endmember extraction,” IEEE Transactions on Geoscience and Remote Sensing, vol. 49, no. 11, pp. 4177–4193, 2011.
  • [16] M. D. Craig, “Minimum-volume transforms for remotely sensed data,” IEEE Transactions on Geoscience and Remote Sensing, vol. 32, no. 3, pp. 542–552, 1994.
  • [17] C. Zhao, G. Zhao, and X. Jia, “Hyperspectral image unmixing based on fast kernel archetypal analysis,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 10, no. 1, pp. 331–346, 2017.
  • [18] G. Zhao, C. Zhao, and X. Jia, “Multilayer unmixing for hyperspectral imagery with fast kernel archetypal analysis,” IEEE Geoscience and Remote Sensing Letters, vol. 13, no. 10, pp. 1532–1536, 2016.
  • [19]

    M. Mørup and L. K. Hansen, “Archetypal analysis for machine learning and data mining,”

    Neurocomputing, vol. 80, pp. 54–63, 2012.
  • [20] A. Cutler and L. Breiman, “Archetypal analysis,” Technometrics, vol. 36, pp. 338–347, 1994.
  • [21] S. Seth and M. J. Eugster, “Probabilistic archetypal analysis,” Machine Learning, vol. 102, no. 1, pp. 85–113, 2016.
  • [22] C. Römer, M. Wahabzada, A. Ballvora, F. Pinto, M. Rossini, C. Panigada, J. Behmann, J. Léon, C. Thurau, C. Bauckhage et al., “Early drought stress detection in cereals: Simplex volume maximisation for hyperspectral image analysis,” Functional Plant Biology, vol. 39, no. 11, pp. 878–890, 2012.
  • [23] C. Seiler and K. Wohlrabe, “Archetypal scientists,” Journal of Informetrics, vol. 7, no. 2, pp. 345–356, 2013.
  • [24] R. Heylen, P. Scheunders, A. Rangarajan, and P. Gader, “Nonlinear unmixing by using different metrics in a linear unmixing chain,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 8, no. 6, pp. 2655–2664, 2015.
  • [25] R. Wang, H. Li, W. Liao, and X. Huang, “Endmember initialization method for hyperspectral data unmixing,” Journal of Applied Remote Sensing, vol. 10, no. 4, p. 042009, 2016. [Online]. Available:
  • [26] M. Zortea and A. Plaza, “A quantitative and comparative analysis of different implementations of n-findr: A fast endmember extraction algorithm,” IEEE Geoscience and Remote Sensing Letters, vol. 6, no. 4, pp. 787–791, 2009.
  • [27] A. Okujeni, S. van der Linden, and P. Hostert, “Berlin-urban-gradient dataset 2009 - an enmap preparatory flight campaign (datasets),” GFZ Data Services, 2016.
  • [28] C. Thurau, K. Kersting, and C. Bauckhage, “Yes we can: simplex volume maximization for descriptive web-scale matrix factorization,” in Proc. International Conference on Information and Knowledge Management.   ACM, 2010, pp. 1785–1788.
  • [29] R. Heylen, D. Burazerović, and P. Scheunders, “Non-linear spectral unmixing by geodesic simplex volume maximization,” Selected Topics in Signal Processing, IEEE Journal of, vol. 5, no. 3, pp. 534–542, 2011.
  • [30] P. J. Green, “Reversible jump markov chain monte carlo computation and bayesian model determination,” Biometrika, vol. 82, no. 4, pp. 711–732, 1995.
  • [31] Z. Zhang, Y. Xu, J. Yang, X. Li, and D. Zhang, “A survey of sparse representation: algorithms and applications,” IEEE access, vol. 3, pp. 490–530, 2015.
  • [32] T. Cocks, R. Jenssen, A. Stewart, I. Wilson, and T. Shields, “The hymap airborne hyperspectral sensor: the system, calibration and performance,” in Proceedings of the 1st EARSeL workshop on Imaging Spectroscopy.   EARSeL, 1998, pp. 37–42.
  • [33] R. Richter and D. Schläpfer, “Geo-atmospheric processing of airborne imaging spectrometry data. part 2: atmospheric/topographic correction,” International Journal of Remote Sensing, vol. 23, no. 13, pp. 2631–2649, 2002.
  • [34] K. Segl, L. Guanter, C. Rogass, T. Kuester, S. Roessner, H. Kaufmann, B. Sang, V. Mogulsky, and S. Hofer, “Eetes—the enmap end-to-end simulation tool,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 5, no. 2, pp. 522–530, 2012.
  • [35] U. Heiden, K. Segl, S. Roessner, and H. Kaufmann, “Determination of robust spectral features for identification of urban surface materials in hyperspectral remote sensing data,” Remote Sensing of Environment, vol. 111, no. 4, pp. 537–552, 2007.
  • [36] C. J. Geyer and J. Møller, “Simulation Procedures and Likelihood Inference for Spatial Point Processes,” Scandinavian Journal of Statistics, vol. 21, no. 4, pp. pp. 359–373, 1994.
  • [37] D. Bulatov, S. Wenzel, G. Häufel, and J. Meidow, “Chain-wise generalization of road networks using model selection,” ISPRS Annals of Photogrammetry, Remote Sensing and Spatial Information Sciences, vol. IV-1/W1, pp. 59–66, 2017. [Online]. Available:
  • [38] M. M. Breunig, H.-P. Kriegel, R. T. Ng, and J. Sander, “Lof: identifying density-based local outliers,” in ACM sigmod record, vol. 29, no. 2.   ACM, 2000, pp. 93–104.