International Journal of Earth & Environmental Sciences Volume 3 (2018), Article ID 3:IJEES-152, 8 pages
Research Article
Wavelets Based Simulation and Visualization Approach for Unmixing of Hyperspectral Data

Mutaz Mohammad1, Fares M. Howari1,*, Gheorge Acbas1, Yousef Nazzal1 and Fatima AlAydaroos2

1College of Natural and Health Sciences, Zayed University, P.O. Box 144534 Abu Dhabi, UAE
2UAE Space Agency, P.O. Box 7133 Abu Dhabi, UAE
Dr. Fares M. Howari, College of Natural and Health Sciences, Zayed University, P.O. Box 144534, Abu Dhabi, UAE; E-mail:
01 June 2018; 28 June 2018; 30 June 2018
Mohammad M, Howari FM, Acbas G, Nazzal Y, AlAydaroos F, et al. (2018) Wavelets Based Simulation and Visualization Approach for Unmixing of Hyperspectral Data. Int J Earth Environ Sci 3: 152. doi:
This study was supported by a Grant-in-Aid for Scientific Research (C) from the Japan Society for the Promotion of Science (grant number 20592601).


This paper discusses the use of wavelet based simulation and visualization methods to unmix multi-component, intimate, and particulate mixtures. We describe a new approach to simulate the hyperspectral data based on an interpolating wavelet estimation. We choose a very constructive family of a wavelet basis, Duabechies wavelets, of order one and two to analyze the data. We used the discrete wavelet transform (DWT), which the basic tool needed for studying time series via wavelets and plays a role analogue in the theory of spectral analysis. The present paper present parameters and algorithms that successfully simulate the binary system of Na2SO4-MgSO4. This system was selected because of it is importance to us on analogue for Mars salt and salt regolith.

1. Introduction

Unmixing hyperspectral signals and prole is still a challenge faces environmental, remote sensing, and geoscience communities. Hyperspectral data provide wealth of information about the physicochemical conditions of earth and extraterrestrial targets. However, the purity level of the spectra is a major challenge. Very often the hyperspectral proles contain several end members or components. In theory if a mixed spectra can be simulated using a mathematical approach then an inversion technique would reveal the end members. This approach was tested to simulate mixed spectra of different ratios of this binary system of Na2SO4-MgSO4. On earth these chemicals salts of high enough concentrations can be considered economic minerals and salts that can be used in several industrial applications. For example, Na2SO4 is used in manufacturing wood pulp, glass industry, thermal storage, drying agent. Whereas magnisum sulfate is used commonly as epsom salt for medical purposes, and used in agriculture to increase sulfer or magnesium concentrations in soil, and it is used as brewing salt in beer among other uses. On Mars these salts have been investigated as an analogue to explorer the presence of microbes that possibly inhabit cold MgSO4 rich brines in Mars. This binary system is also an adequate analogue of the martian salts and its salty regolith. For more details about the endmember extraction and unmixing, we refere to [8-14,18-21].

An example from spectral prole extracted from AVRIS image of White Sand in New Mexico, USA is shown in Figure 1. However, the data used in simulation here were obtained using the HR-1024i from the Spectra Vista Corp. (SVC) which is the latest model from their next generation of high performance single-beam eld spectroradiometer measuring over the visible to short-wave infrared wavelength range (350-2500nm).

figure 1
Figure 1: Spectral proles matching as extracted from AVRIS image.

Reflectance spectroscopy and hyperspectral imaging analysis has picked up because it is reliable, fast, and less expensive and not intrusive [7]. However, spectral patterns for mixed materials can't be visually understood. Spectral pretreatment techniques must be applied to smooth spectral graph such as data normalization, continuum removal among other methods. We hypothesis that if the mixed spectra can be simulated mathematically, then we can obtain the end members that make up the mixed spectra. Toward this end, mathematical approaches or models that describe the reectance process in terms of several variables that control light reflection have been used extensively by Hapke in [6]; Pieters and Mustard in [1,15]; Robertson et al. in [16]; and Grumpe et al. in [5]. In the present study, we use wavelet approach to simulate and invert the mixed spectra. Wavelet methods are simple and computationally effective, and can be implemented in real-time. It is proven that wavelet reconstruction is an efficient way to represent functions, operators and big data set due to the capability of wavelet coecients to characterize image/signal discontinuities (i.e., noise) at different scales. In fact discrete wavelet transform (DWT) can be used in various applications, such as data simulation, image compression and coding. DWT refers to wavelet transforms for which the wavelets are discretely sampled. Wavelets provide a spatial frequency reconstruction, very useful in smoothing problems, in particular in density and regression approximation, having excellent statistical properties in data smoothing. They oer frequency and location time representation of data allowing adaptive ltering, estimation and smoothing. Two of the main advantages of wavelet representation is to provide an interpretation of the spectra and minerals related to the given data, using few number of high pass lters (coeffcients) and to provide a more concise representation because it minimizes the amount of redundancy in the coefficients and used to remove sparse noise from the signal.

Information of a hyperspectral image is heavily related to the shape of reectance spectra, which is recovered and represented in the magnitudes of its wavelet coefficients. The wavelet transform is an effective tool in many image processing applications due to the capability of wavelet coeffcients to characterize image/signal discontinuities (i.e., noise) at different scales. The DWT provides a more concise representation because it minimizes the amount of redundancy in the coefficients and used to remove sparse noise from the signal.

The simulation is done by performing three methods, Lagrange interpolation, and Daubechies and Haar wavelet threshold denoising methods, to simulate the given data. The original data can be recovered rst from Lagrange coffiecients and secondly by the (Daubechies and Haar) wavelets coffiecients through so called discrete threshold wavelet lter banks. In order to facilitate the simulation analysis and processing, and as a sample, we considered the reflectances of the wave lenght range 250-2500 nm of the binary system of Na2SO4-MgSO4, then its Lagrange and wavelet threshold transform/interpolation should be discretized. The advantages of the interpolating curves using the methods above is that the entire data can be inferred as if we observe a dense enough set of points and only small amount of coffiecients are needed to achieve accurate approximation.

2. Metheodolgy

In this study, laboratory experiments under controlled conditions have been carried out to prepare pure Na2SO4 and MgSO4 crusts and their mixtures. Analytical grade compounds of Na2SO4 and MgSO4 were considered in this study. The salts were dissolved in water in 1000-ml volumetric asks. Total quantity of 200 ml of each solution were removed and placed in glass Petri dishes. The water was removed from the Petri dishes by evaporation at 40° C in an electric oven for 12 hours, then the samples were removed from the oven. Immediately following the removal from the oven, the samples were placed in a desicator, after which the reectance readings were made. A HR-1024i spectroradiometer was used to obtain the reflectance readings. Most of salts encountered in soils are a mixture of two or more type of salts. Sets of chemically mixed salt samples were prepared from which the pure salts were made using the same steps as described above. The mixing proportions are in 1 : 0, 0:75 : 0:25, 0:5 : 0:5, 0:25 : 0:75, and 0 : 1 weight mixing ratios. The spectra of the mixed samples were compared with the spectra of the pure samples solution fraction. All the pure and mixed samples were examined under petrographic and binocular microscopes for crystal size and morphology observations. Another runs were done by placing soil in the Petri dishes allowing the salt crus to grow on the top. In the ternary system consisting of sodium sulfate and magnesium sulfate, the following solid phases appear in the temperature range from -10° to 110° C:

  1. Ice
  2. Na2SO4∙10H2O, glauber's salt
  3. Na2SO4, thenardite
  4. MgSO4∙12H2O, magnesium sulfate dodecahydrate
  5. MgSO4∙7H2O, epsom salt
  6. MgSO4∙6H2O, hexahydrite
  7. MgSO4∙H2O, kieserite
  8. Na2SO4∙MgSO4∙4H2O, bloedite
  9. Na2SO4∙MgSO4∙2.5H2O, löweite
  10. Na2SO4∙MgSO4, vanthoffite

The temperatures and concentration ranges at which these solids appear are shown in the Figure 2. The phase diagrams shown on these pages are calculated with the Extended UNIQUAC thermodynamic model ( The equilibrium lines and the experimental data in the diagram represent compositions and temperatures at which two solid phases are in equilibrium with the same solution. The water content of the solutions are not shown. It can be thought of as a third dimension in the diagram. The green lines in the gure above are tie lines indicating phases and compositions in equilibrium with each other.

figure 2
Figure 2: Phase diagram for The Na2SO4-MgSO4-H2O system.

The study will apply Lagrange interpolation and wavelet coffiecient analysis. The Lagrange interpolating polynomial is the polynomial P(x) of degree ≤ n-1 that passes through the n points (xk, f(xk)), for k=1,...,n and is dened by

P(x)= k=1 n P k (x)


P k (x)= y k Π j=1,kj n x x j x k x j

When constructing interpolating polynomials, there is a tradeoff between having a better fit and having a smooth well-behaved fitting function. The more data points that are used in the interpolation, the higher the degree of the resulting polynomial, and therefore, a highdegree interpolation may be a poor predictor of the function between points, although the accuracy at the data points will be "perfect".

Daubechies was first to construct compactly supported orthogonal wavelets with a preassigned degree of smoothness [2]. She intended to construct a wavelet with N vanishing moments and supported in [0; 2N-1]; where a wavelet function ψ is said to have N (≥ 2) vanishing moments if,

x n ψ(x)dx =0, n=0,1,...,N1

When N = 1, it is Haar wavelet and has only one vanishing moment. In fact, all systems built by using the unitary extension principle (UEP) of Ron and Shen [17] have only one vanishing moment.

Daubechies looked for hk's in a dilation equation,

ϕ= 2 k h k ϕ(2xk)

such that the orthonormal condition ϕ ^ (2pπ)= δ 0,p , where p = 0,1,..., N-1, is satisfied, and hk are as follows:

h 0 = υ(υ1) 2 ( υ 2 +1) , h 1 = (υ1) 2 ( υ 2 +1) , h 2 = (υ+1) 2 ( υ 2 +1) , h 3 = (υ+1) 2 ( υ 2 +1)

By choosing υ= ±1 3 , a scaling function ϕ has support [0, 2] :

Daubechies wavelet does not have a closed form, but instead, can be obtained recursively by

ϕ n (x)= k 2 h k ϕ n1 (2xk)

ϕ 0 (x)= X [0,1) (x)

It is known that the smoothness of the wavelets increases with N. For an application in numerical analysis, Coifman asked Daubechies to construct a family of wavelets ψ having N vanishing moments, minimum size support and


in such a way we have a smooth orthogonal system [2]. Daubechies wavelet has vanishing moment for wavelet function ψ only. However, Daubechies designed, at that request, a wavelet (Coiets) that has vanishing moment for both wavelet and scaling functions ψ and ϕ. The wavelet is near symmetric such that wavelet function ψ has 2N vanishing moments and the scaling function has 2N-1 vanishing moments. The two functions have a support width of length 6N- 1. Also there is no closed form for Coiets and it can be obtained recursively.

2.1 Denition

[3] A compactly supported function ϕ L 2 () is said to be renable if

ϕ(x)=2 k h 0 [k]ϕ(2xk)

for some nite supported sequence h 0 [k] l 2 () . The sequence h0 is called the low pass filter of ϕ.

Let Ψ= { Ψ l } l=1 r L 2 () be of the form

Ψ l =2 k h 1 [k]ϕ(2k)

where { h 1 [k],k} l=1 r is a nitely supported sequence and is called a high pass filter of the system. Dene the wavelet system

X(Ψ)={ Ψ l,j,k :1lr;j,k}

where Ψ l,j,k (x)= D j T k Ψ l . The general setup is to construct a wavelet for L 2 () and of the form of X(Ψ) [4].

The Haar function is dened by

ϕ= 1 [0,1);

V j ={f L 2 ():f  is constant on [2 j k, 2 j (k+1])], k}

Note that the Haar wavelet can be written as,

Ψ(x)= X [0,1/2) (x) X [1/2,1) (x)

= X [0,1) (2x) X [0,1) (2x1)

= 1 2 (D X [0,1) (x)D T 1 X [0,1) (x))

Already in 1910 it was proved by Haar that the functions Ψ j,k ,j,k constitute an orthonormal basis for L 2 () . The Daubechies wavelets functions for different order are given in Figure 3.

figure 3
Figure 3: Daubechies scaling functions of order 1,2 and 3, respectively.

3. Results and Discussion

3.1 Lagrange and Wavelet Simulations

To illustrate the simulation, assume that X(Ψ) forms a wavelet system for L 2 () . For any signal/data f L 2 (),f can be expanded as follows,

f= l=1 r j k f, Ψ l,j,k Ψ l,j,k

Our simulation is based on signal classication approach using discrete wavelet transform (DWT). A DWT-based linear (un) mixing system is designed specially for estimation. Using the wavelet transform, the original signal is represented by a set of wavelet transform coffiecients, and features are extracted from these coffiecients. Wavelet theory provides a good technique for signal approximation using scaled basis functions. The inverted data of the hyperspectral imaging illustrated in Figure 4. Note that the XY-axes ranges presented in Figures 4-11 are not refer to the actual wavelength and reflectance ranges of this study. The ranges presented here came from the ListLinePlot command using Mathematica software that plots a line through the points {1, y1}, {2, y2},..., {647, y647} where the set {y1, y2,..., y647} presenting the spectral image data of the binary system of Na2SO4-MgSO4.

figure 4
Figure 4: Hyperspectral imaging (un)mixing data.
figure 5
Figure 5: Simulation using Lagrange interpolation.
figure 6
Figure 6: A closer view of the simulation using Lagrange interpolation.
figure 7
Figure 7: Simulations using Haar and Daubechies wavelets thresholds, respectively.
figure 8
Figure 8: A closer view for the simulations using Haar and Daubechies wavelets thresholds, respectively.
figure 9
Figure 9: Haar wavelet simulation of the whole spectra.
figure 10
Figure 10: Splitting the the binary system of Na2SO4-MgSO4 using the vector projection method.
figure 11
Figure 11: Splitting the whole spectra using the vector projection method.

The first approach utilized the Lagrange Interpolation Technique. The simulation using Lagrange interpolation is illustrated in Figure 5. A closer view of the simulation using this approach is depicted in Figure 6.

The second approach is by using wavelets interpolation, the first will be done using Haar wavelets and the second by Daubechies wavelets. The simulations using these functions are depicted in Figure 7 and 8.

3.2 Splitting the spectra using vector projection

For this section, we represent the hyperspectral data as a sum of two functions that makes it. One way doing that is to use the vector projection method. The vector projection method of a vector a on (or onto) a nonzero vector b is the orthogonal projection of a onto a straight line parallel to b . It is a vector parallel to b , defined as

pro j b a = a b b b b

To illustrate the idea, we choose a specic vector b such that a = c + d

where c = proj b a and d = a proj b a .

Figure 10 shows the graphs of these splitted simulations.

4. Conclusion

The study used two approaches 1) semantic one by tracing the diagnostic spectral features such as the location, shape and depth of the band for certain endmembers and through study the association of the diagnostic feature between the mixed spectra and their endmembers (semantic or diagnostic approach), 2) mathematically through wavelet approach that simulate and unmix the reflectance spectra. To this end, a one-dimensional threshold method for hyperspectral images has been proposed based on Lagrange interpolation and discreet wavelet threshold denoising method. The performance of Lagrange interpolation method and the wavelet approach based on Haar and Daubechies threshold methods to simulate Na2SO4-MgSO4 binary system yielded satisfactory results. Vector projection method was used successfully to split the hyperspectral data as a sum of two data functions that simulate it.

5. List of Symbols

P(x) A polynomial of degree n
The set of real numbers
The set of integers, 0,±1,±2,...
h0[k] Sequence of real number where k
k=1 n x k The summation of the values of xk where k
k x k Innite summation
Ψ The wavelet function
f(x)dx The integral operator of the function f
j=1 n x k Finite product
L 2 () The space of functions f: for which f 2 <
Taf(x) The translation operator by a of the function f
Df(x) The dilation operator of the function f dend by 2 f(2x)
Ψ l,j,k Ψ l,j,k (x)= D j T k Ψ l (x)= 2 j/2 Ψ( 2 j xk)
a∙b The dot product of a and b
[a,b] The set of real numbers between a and b
X[0,1)(x) The characteristic function of the set [0, 1]
f, Ψ l,g,k The inner product between f and Ψ l,g,k
X(Ψ) The wavelet system

Competing Interests

The authors declare that they have no competing interests.


We are very grateful to the anonymous referees valuable comments and suggestions.


  1. Cheek LC, Pieters CM (2014) Re ectance spectroscopy of plagioclasedominated mineral mixtures: Implications for characterizing lunar anorthosites remotely. American Mineralogist 99: 1871-1892. View
  2. Daubechies I (1988) Orthonormal bases of compactly supported wavelets. Comm Pure Appl Math 41: 909-996. View
  3. Daubechies I (1992) Ten Lectures on Wavelets. SIAM, Philadelphia, PA.
  4. Daubechies I, Han B, Ron A, Shen Z (2003) Framelets: MRA-based constructions of wavelet frames. Appl Comput Harmon 14: 1-46. View
  5. Grumpe A, Mengewein N, Rommel D, Mall U, Wöhler C, et al. (2018) Interpreting spectral unmixing coe cients: From spectral weights to mass fractions. Icarus 299: 1-14. View
  6. Hapke B (1981) Bidirectional re ectance spectroscopy: 1. Theory J of Geop Resh 86: 3039-3054.
  7. Howari FM, Banat KM (2001) Assessment of Fe, Zn, Cd, Hg, and Pb in the Jordan and Yarmouk River Sediments in Relation to Their Physicochemical Properties and Sequential Extraction Characterization. Water, Air, Soil Pollution 132: 43-59. View
  8. Howari FM, Goodell PC (2009) Characterization of salt-crust buildup and soil salinization in the United Arab Emirates by means of field and remote sensing techniques. Remote sensing of soil salinization: Impact on land management.
  9. Howari FM (2009) Remote Sensing of Soil Salinization in Soil salinity and salinization hazard/ A Zinck. G Metternicht Ed.
  10. Howari FM, BR Jordan, N Bouhouche, Wyllie-Echeverria S (2009) Field and remote-sensing assessment of mangrove forests and seagrass beds in the northwestern part of the United Arab Emirates. J of Coa Res 25: 48-56. View
  11. Howari FM (2006) Spectral analyses of sabkha sediments with implications for remote sensing on Mars. Int J of Astro 5: 47-56. View
  12. Howari FM (2004) Chemical and Environmental Implications of Visible and Near-Infrared Spectral Features of Salt Crusts Formed from Di erent Brines. Annali di chimica 94: 315-323. View
  13. Howari FM (2003) Comparison of spectral matching algorithms for identifying natural salt crusts. J of app spec 70: 782-787. View
  14. Howari FM, Goodell PC, Miyamoto S (2002) Spectral properties of salt crusts formed on saline soils. J of Envir Qual 31: 1453-1461. View
  15. Mustard JF, Pieters CM (1989) Photometric Phase Functions of Common Geologic Minerals and Applications to Quantitative Analysis of Minerals Mixture Reflectance Spectra. J of Geop Res 94: 13,619-13,634. View
  16. Robertson KM, Milliken RE, Li S (2016) Estimating mineral abundances of clay and gypsum mixtures using radiative transfer models applied to visible-near infrared reflectance spectra. Icarus 277: 171-186. View
  17. Ron A, Shen Z (1997) Affne systems in L2(Rd): the analysis of the analysis operators. J Funct Anal 148: 408-447. View
  18. A Sadiq, Howari F (2006) Remote sensing and spectral characteristics of desert sand from Qatar Peninsula, Arabian/Persian Gulf. Remote Sensing 1: 915-933. View
  19. Stagakis S, Vanikiotis T, Sykioti O (2016) Estimating forest species abundance through linear unmixing of CHRIS/PROBA imagery. ISPRS J of Photo and Remo Sen 119: 79-89. View
  20. Sun W, Yang G, Wu K, Li W, Zhang D, et al. (2017) Pure endmember extraction using robust kernel archetypoid analysis for hyperspectral imagery. ISPRS J of Photo Remo Sen 131: 147-159. View
  21. Zhang C, Qin Q, Zhang T, Sun Y, Chen C, et al. (2017) Endmember extraction from hyperspectral image based on discrete re y algorithm (EEDFA). ISPRS J of Photo Remo Sen 126: 108-119. View