Montmorillonite Estimation in Clay–Quartz–Calcite Samples from Laboratory SWIR Imaging Spectroscopy: A Comparative Study of Spectral Preprocessings and Unmixing Methods

Abstract

Clay minerals play an important role in shrinking–swelling of soils and off–road vehicle mobility mainly due to the presence of smectites including montmorillonites. Since soils are composed of different minerals intimately mixed, an accurate estimation of its abundance is challenging. Imaging spectroscopy in the short wave infrared spectral region (SWIR) combined with unmixing methods is a good candidate to estimate clay mineral abundance. However, the performance of unmixing methods is mineral-dependent and may be enhanced by using appropriate spectral preprocessings. The objective of this paper is to carry out a comparative study in order to determine the best couple spectral preprocessing/unmixing method to quantify montmorillonite in intimate mixtures with clays, such as montmorillonite, kaolinite and illite, and no-clay minerals, such as calcite and quartz. To this end, a spectral database is built with laboratory hyperspectral imagery from 51 dry pure mineral samples and intimate mineral mixtures of controlled abundances. Six spectral preprocessings, standard normal variate (SNV), continuum removal (CR), continuous wavelet transform (CWT), Hapke model, first derivative (1st SGD) and pseudo–absorbance (Log(1/R)), are applied and compared with reflectance spectra. Two linear unmixing methods, fully constrained least square method (FCLS) and multiple endmember spectral mixture analysis (MESMA), and two non-linear unmixing methods, generalized bilinear method (GBM) and multi-linear model (MLM), are compared. Global results showed that the benefit of spectral preprocessings occurs when spectral absorption features of minerals overlap for SNV, CR, CWT and 1st SGD, whereas the use of reflectance spectra performs the best when no overlap is present. With one mineral having no spectral feature (quartz), montmorillonite abundance estimation is difficult and gives RMSE higher than 50%. For the other mixtures, performances of linear and non-linear unmixing methods are similar. Consequently, the recommended couple spectral preprocessing/unmixing method based on the trade-off between its simplicity and performance is 1st SGD/FCLS for clay binary and ternary mixtures (RMSE of 9.2% for montmorillonite–illite mixtures, 13.9% for montmorillonite–kaolinite mixtures and 10.8% for montmorillonite–illite–kaolinite mixtures) and reflectance/FCLS for binary mixtures with calcite (RMSE of 8.8% for montmorillonite–calcite mixtures). These performances open the way to improve the classification of expansive soils.

1. Introduction

Soil shrinking and swelling due to expansive clays have consequences on urban planning for decision makers since they cause billions of dollars in building damages every year [1,2,3]. They are also responsible for off-road vehicle mobility, the latter being impacted by the sinking and stickiness of wet soils [4,5,6]. Aluminum-rich smectite (montmorillonite) has been demonstrated to be the clay mineral having the most expansive potential hazard [7]. Thus to monitor soil shrinking and swelling, the characterization of clay mineralogy is needed, both in terms of detection and quantification. In mid-latitude regions, smectite (whose montmorillonite represents the most common variety), illite and kaolinite are the main clay minerals found in soils affected by swelling risk [8]. For soil swell–shrinking application, Dufréchou et al. and Chassagneux et al. [9,10] defined four classes linking shrink–swell potential and montmorillonite content: low swelling potential (<10%), moderate swelling potential (between 10% and 50%), high swelling potential (between 50% and 70%) and very high swelling potential (>70%). Many well-used techniques, either qualitative or quantitative, can assess clay minerals characterization, but they only provide little and local information due to a low number of soil samplings along a coarse grid such as: ground geotechnical engineering techniques measuring soil swelling potential [11] and X-ray diffraction (XRD) for mineral identification/quantification, and infrastructure damage reports for expansive clay qualitative assessment. However, the use of these methods mainly remains expensive and time-consuming.

An alternative to overcome these limitations is the use of hyperspectral imaging spectroscopy, which is able to discriminate clay minerals using their specific spectral absorption features in the short wave infrared region (SWIR) 2100–2500 nm [12,13,14]. However, clay minerals are rarely pure in soils and are usually intimately mixed with other minerals, water and organic matter. Spectral preprocessings are used to overcome intraclass spectral variability, remove some non-linear illumination effects and enhance shallow absorption features [15,16,17,18,19,20]. For example, the identification of soil minerals using spectral similarity distance increases up to 18% using continuous wavelet transform [20]. Using spectral preprocessings could be a first preliminary step to further discrimination/quantification methods [16,20,21].

Two classes of methods exist in order to map mineral clays in soils: the first is based on spectral features and the second on unmixing. The first class includes the use of spectral indices (Dufréchou et al. [12] obtained a root mean square error (RMSE) of 15.4% for montmorillonite, 25.2% for kaolinite and 29.8% for illite), and linear regression methods using clay spectral feature characteristics (e.g., depth, width and position close to 2200 nm; [22,23,24,25]), and decision tree methods (Mulder et al. [26] obtained a RMSE of 9% for clay–calcite intimate mixtures). Expert systems such as Tetracorder use spectral features in order to detect and discriminate minerals, such as clays [27], for example to map clay alteration of Hawaii volcanoes [28]. These characteristics can also be automatically extracted using partial least square regression (Viscarra Rossel et al. [29,30] obtained a RMSE of 3% for illite, kaolinite and smectite in clay–quartz–organic mixtures). Neural networks have been tested for mineral sub-pixel classification [31], on a Hyperion dataset. Compared to ground reference, fraction prediction provides a correlation coefficient r²; value of 0.62 for illite/montmorillonite class and 0.6 for kaolinite class.

The second class of methods using unmixing relies on the knowledge of pure mineral spectra or endmembers [32,33,34,35]. Linear unmixing methods only consider materials at the surface that are contiguously distributed (named in the following “aeral mixture”, [36] Figure 1a) while non linear unmixing methods consider mixed pixels several soil elements of volumetric distribution (also called “intimate mixture”, [32,36] Figure 1b). On one side, the most used linear unmixing methods considering only one spectrum per endmember class are the fully constrained least square (FCLS; [37]) and mixture tuned matched filtering (MTMF; [38]). On contrary, several algorithms account for the endmember spectral intraclass variability such as multiple endmember spectral mixture analysis (MESMA; [39]), the spectral assistant, (TSA; [40]) and the unmixing within a multi-task Gaussian process framework (UMTGP; [41]). For example, Chabrillat and Bedini [42,43] used MTMF and MESMA to discriminate illite/smectite from kaolinite in soils with airborne data. On the other side, Heylen [44] gave a review of performance of nonlinear unmixing methods, including the Hapke model [45], generalized bilinear model (GBM; [46]) and multi-linear model (MLM; [34]). For example, [34] showed with laboratory spectroscopy datasets of quartz–alunite intimate mixtures that the Hapke model delivers the best result (bias of 1–2% between estimated and measured abundances), then MLM (bias of 10–20%) and at last both GBM and FCLS have a similar performance (bias 30–40%). Robertson et al. [47] compared a version of the Shkuratov model [48] with the Hapke model to quantify laboratory mineral mixtures of montmorillonite–gypsum. In these samples, minerals fractions are estimated with an error less than 10% for both models. Unmixing methods can also be embedded in expert systems such as Tetracorder in order to increase the accuracy of mineral discrimination [27].

The advantage of using unmixing methods is to be not site-dependent and they do not require a learning stage with a representative training dataset, which is the case with the first class of methods. The difficulty to choose the most appropriate unmixing method is that their performance differs from the mineral mixture composition of soils. Few studies have used unmixing methods to estimate clay minerals fraction, but unmixing is more used for mineral discrimination. Moreover, very few studies have combined spectral preprocessings with unmixing methods, but they demonstrated the ability to decrease mineral fraction estimation errors using spectral preprocessing [16,20,21].

Then, by using airborne and satellite hyperspectral imaging spectroscopy, the sensor spatial resolution and the possible presence of vegetation overlaid on soil within a pixel degrade the performances of mineral estimation [42,43,49]. A first step commonly practiced is to study spectroscopy data under laboratory conditions to avoid environmental factors such as atmospheric and soil water content and soil/vegetation mixtures [12,26,29,50,51,52].

The main objective of this study is to compare the performance of both linear and non-linear unmixing methods combined with spectral preprocessings to estimate montmorillonite abundance in mineral mixtures. For this purpose, this pioneer work is based on spectroscopic measurements of intimate mixtures composed of different mineral types (clays, calcite and quartz) and controlled abundances, manually generated in the laboratory with dry conditions and further used as proxies for soil samples. Performances of the spectral preprocessings and unmixing methods will be assessed from the spectral database deriving from them. Another scenario considering synthetic mineral aeral mixtures computed from spectral measurements of pure minerals will be also analyzed in order to compare the performances of spectral preprocessings and unmixing methods between aeral and intimate mixtures. Thus, the spectral database of intimate and aeral mixtures and the laboratory imaging spectroscopy setup are presented in Section 2. The methodology including the spectral preprocessing and selected unmixing methods are described in Section 3. Results are presented in Section 4 and discussed in Section 5, with the conclusions exposed in Section 6.

2. Material

2.1. Laboratory Imaging Spectroscopy Setup

The laboratory setup (Figure 2) was similar to the one used in [50] under dry conditions. Two halogen lamps light the sample with an incidence zenith angle of 35°. The reflected signal was acquired by a hyperspectral camera HySpex SWIR–320me [53] (NorskElektroOptikk) whose characteristics are summarized in

Table 1. Characteristics of the HySpex hyperspectral camera.

	Hyspex SWIR–320m–e Camera
Spectral range (nm)	1000–2500
Spectral resolution (nm)	6
Number of spectral bands	256
Field of view (°)	13.5
Number of pixels across track	320

. Each dry mineral mixture was placed in a plastic container lying on a tray that moved horizontally during the acquisition process to obtain a hyperspectral image Rad_sample (Figure 2). The final image was the average of 10 acquisitions, reducing the noise by $\sqrt{10}$. Further, an image, Rad_spec, of a standard white Spectralon® panel was recorded in order to convert the acquired sample images expressed in radiance unit to reflectance unit such as:

$$\rho {\left(\lambda ,i,j\right)}_{sample}=\frac{Rad{\left(\lambda ,i,j\right)}_{sample}}{Rad{\left(\lambda ,i,j\right)}_{spec}}\times \rho {\left(\lambda \right)}_{spec}$$

(1)

With Rad(λ, i, j) the radiance of image pixel (i, j) with i the line number, j the column number and λ the wavelength, ρ(λ)_spec the Spectralon^® reflectance and ρ(λ,i,j)_sample the sample reflectance.

An area of 150 pixels × 200 pixels (11.25 × 15.0 cm²) was selected over each image to avoid possible stray light from the container borders (Figure 2b). Saturated pixels were removed, corresponding to 2% in the average for each image. For each mineral mixture, a mean reflectance over the area of interest of the image was computed and stored in a spectral database.

2.2. Spectral Database

The spectral database was built from five pure clay mineral samples (illite, kaolinite, montmorillonite, calcite and quartz) and two mixture types were considered, aeral and intimate, for a total of 46 mixtures. In one hand, the spectral reflectance of aeral mixtures was computed from the linearly mixed combination of pure mineral samples pondered by their respective abundance. On the other hand, the spectral reflectance of intimate mixtures was measured from the facilities presented in Section 2.1. This last scenario is a proxy to consider soil mixtures, without any organic matter and soil water content affecting spectral signatures [54]. Investigating unmixing performance with aeral mixtures enables to (i) highlight the benefit or not of spectral preprocessings over reflectance spectra, and (ii) compare the performance of linear unmixing methods independently of non-linear effects occurring in intimate mixtures. For this purpose, the mineral abundance should be as close as possible between aeral and intimate mixtures.

Each sample was labeled as a suite of capital initials of the mineral’s name followed by its abundance in weight percentage (wt %), and the same for quartz and calcite. For example, I20K80 refers to illite and kaolinite with an abundance of 20 wt % and 80 wt %, respectively. M20S70C10 corresponds to a mixture with a 20 wt % abundance of montmorillonite, 70 wt % of quartz and 10 wt % of calcite, respectively.

2.2.1. Pure Clay Mineral Samples

The following clay minerals are: illite from Saint Paulien (France) commercialized by Argiles du Velay, kaolinite from Clérac (France) commercialized by Imerys, montmorillonite from Villanova Tulo (Sardinia, Italy) commercialized by Argiles du Bassin Méditerranéen (ABM). Clays are provided as a powder with a grain size less than 80 µm. Quartz is collected from a sand mine in Fontainebleau (France) commercialized by Sifraco and calcite comes from Orgon (France). Quartz has a grain size less than 300 µm, and calcite less than 70 µm. Collected minerals were oven-dried and have been stored in dry conditions before use.

Clay spectral features are present in the SWIR domain, at 1400 nm, 1900 nm and more specifically, between 2100 and 2500 nm [12,13,14] (Figure 3). As atmospheric water vapor content has low transmittance at 1400 and 1900 nm [14,55], only the reduced spectral domain of 2100 – 2425 nm will be further considered in this study. Clays exhibit several spectral features, which result from the presence of vibrational hydroxyl processes:

• The absorption feature around 2200 nm is due to the Al–OH vibrational mode. Its accurate location depends on the clay type: 2208 nm for illite, 2212 nm for montmorillonite and 2206 nm for kaolinite. Kaolinite also has a double absorption feature (2160 nm and 2206 nm), which is leftward asymmetric. The typical absorption bandwidth is around 100 nm whatever the clay type.
• The absorption features due to OH-stretching bands combined with lattice vibrations at approximately 2360 nm is shallow for illite and sharp for kaolinite. Kaolinite has also in addition two more absorption features at 2320 nm and 2380 nm.

Calcite and quartz were very reflective with a maximum reflectance of 0.94 for calcite and 0.8 for quartz (Figure 3). Calcite had strong absorption bands corresponding to the vibrational modes of ${\mathrm{CO}}_3^{2-}$ions between 2300 and 2350 nm, and a weak absorption feature close to 2140 nm. Quartz remained spectrally flat in this spectral domain.

2.2.2. Synthetic Mineral Aeral Mixtures

Aeral mixtures are synthetically generated using spectral data of pure minerals (Figure 3). In order to reproduce variability of pure mineral spectra, 1000 pixel spectra have been randomly selected in each hyperspectral image of pure minerals. Then, using these spectra, a linear mixing was applied to these spectra for each mineral following the equation below:

$$\rho ={{\displaystyle \sum }}_{i=1}^p\left(f_i{\rho }_i\right) with 0\le f_i\le 1 \mathrm{and} {{\displaystyle \sum }}_{i=1}^p\left(f_i\right)=1,$$

(2)

With $f_i$ the abundance of mineral i and ${\rho }_i$ the spectral reflectance. As the focus of the study is on the discrimination of the montmorillonite, only binary aeral mixtures combined with montmorillonite are considered (Figure 4): illite–montmorillonite, kaolinite–montmorillonite, calcite–montmorillonite and quartz–montmorillonite. Each mineral abundance of a given mixture was measured by weight percentage (wt %). The abundances vary from 20 to 80 wt % by a step of 20%.

2.2.3. Mineral Intimate Mixtures

Twelve binary and 9 ternary mixtures of clay minerals (Figure 4a) and 8 binary mixtures of montmorillonite with either quartz or calcite (Figure 4b) were considered in order to uniformly cover as much as possible the area of the ternary diagrams. Each mixture had a volume of 300 mL (corresponding to a sample thickness of 2 cm in a plastic container of surface 22 × 26 cm²). This abundance was converted from volume to weight percentage in order to be compared with its aeral counterpart. As it can be seen in Figure 4, this measurement protocol did not allow us to necessarily have a superposition in terms of mineral abundance between the intimate and patchwork like mixture.

Then, the minerals were manually and homogeneously mixed in the plastic container. The surface of each sample was leveled to avoid non-linear interactions due to the roughness of the sample. The mixture homogeneity was assessed by deriving the standard deviation of the spectral reflectance measured over the image area of interest (cf. Figure 2, Section 2.1)). The latter ranged between 1% to 2.5% reflectance for each mixture, which proved good uniformity within the sample at the surface.

3. Methods

The methodology consists of the application of spectral preprocessing (Section 3.1), then unmixing methods (Section 3.2), and the accuracy of the retrieved abundance maps was further assessed with evaluation criteria (Section 3.3; Figure 5).

3.1. Spectral Preprocessings

Six spectral preprocessings already used with unmixing methods or commonly used in the soil scientific community were selected (

Table 2. Selected spectral preprocessing methods: reflectance at wavelength λ:ρ_λ, averaged reflectance: $\overline{\rho }$, wavelet-base function: ψ, wavelet scale: S, cosinus of the incident zenithal angle: μ₀, cosinus of the viewing zenithal angle: μ, the half width of the smooth filter: m, smoothed spectra: ${\overline{\rho }}_j$, center of the smoothing window: ${\rho }_j$, coefficient of the $i$th point of the filter of the $q$th order of derivative: $P_i^{\left(q\right)}$.

Preprocessing	Equation	Reference
Pseudo absorbanceLog(1/R)	${\rho }_{\lambda }{}^{\prime }=\log \left(1/{\rho }_{\lambda }\right)$	[19] (Rinnan et al., 2009)
Hapke Model (W)	${\rho }_{\lambda }{}^{\prime }=1-{\left(\frac{\sqrt{{\left({\mu }_0+\mu \right)}^2{\rho }_{\lambda }{}^2+\left(1+4\mu {\mu }_0{\rho }_{\lambda }\right)\left(1-{\rho }_{\lambda }\right)}-\left({\mu }_0+\mu \right){\rho }_{\lambda }}{\left(1+4\mu {\mu }_0{\rho }_{\lambda }\right)}\right)}^2$	[44] (Heylen et al., 2014)
Standard Normal Variate (SNV)	${\rho }_{\lambda }{}^{\prime }=\frac{\left({\rho }_{\lambda }-\overline{\rho }\right)}{\sqrt{\frac{\sum {\left({\rho }_{\lambda }-\overline{\rho }\right)}^2}{n-1}}}$	[59] (Barnes et al., 1989)
Continuum Removal (CR)	${\rho }_{\lambda }{}^{\prime }={\rho }_{\lambda }/Continuu{m}_{\lambda }$	[60] (Clark and Roush, 1984)
Continuous Wavelet Transform (CWT)	${\rho }_{\lambda }{}^{\prime }={\displaystyle {\displaystyle \sum }_{{\mathsf{\lambda }}^{\prime }=0}^{N-1}}{\rho }_{\mathsf{\lambda }}, \psi \left[\frac{({\mathsf{\lambda }}^{\prime }-\mathsf{\lambda })\mathsf{\delta }\mathrm{t}}{S}\right]$	[20] (Rivard et al., 2008)
First Savitzky–Golay Derivative (1st SGD)	${\rho }_{\lambda }{}^{\prime }=\frac{d^q{\overline{\rho }}_j}{d{\lambda }^q}={\displaystyle {\displaystyle \sum }_{i=-m}^m}{P}_i^{\left(q\right)}{\rho }_{j+i}$	[61] (Tsai and Philpot, 1998)

). Among them, a first group conserving the reflectance continuum includes pseudo-absorbance (Log(1/R)), the Hapke model and standard normal variate (SNV), while a second group suppressing it includes continuum removal (CR), continuous wavelet transform (CWT) and first Savitzky–Golay derivative (1st SGD).

The chosen Hapke model considers the soil medium as an isotropic mixture with the same granulometry for all components [45]. The acquisition conditions were fixed to 35° for the zenithal incident angle (i.e., μ₀ =cos(35)) and 0° for the zenithal viewing angle (i.e., μ =cos(0)).

CWT deconvolves the signal into a linear sum of a wavelet-base function (ψ) defined as the second order derivative of Gaussian. Ten wavelets were computed from the original spectrum. Consequently, only the second to the fifth scale wavelets were kept in order to remove the contribution of the noise and continuum as recommended by Feng et al. [56].

The 1st SGD is the result of the convolution of the reflectance spectra by the Savitzky–Golay filter [57]. It was performed with a second order polynomial fit on a window size of seven bands. Higher orders of derivative were not tested since Heydley et al. [58] showed that they did not improve unmixing performance.

3.2. Unmixing Methods

Four unmixing methods were chosen: two linear methods, FCLS (the most commonly used) and MESMA (accounts for intraclass variability), and two non-linear methods, GBM (able to considered first and second order interactions) and MLM (multiple interactions;

Table 3. Unmixing methods. Reflectance: ρ, number of endmembers: p, endmember abundance: ${\mathrm{f}}_{\mathrm{i}}$, endmembers spectra: ${\mathsf{\rho }}_{\mathrm{i}}$ and ${\mathsf{\rho }}_{\mathrm{j}}$, residual term: $\mathsf{\epsilon }$, free additional parameter ${\mathsf{\gamma }}_{\mathrm{ij}}$, probability to undergo multiple interactions: P.

Unmixing	Equation	Reference
Fully Constraint Least Square (FCLS)	$\rho ={\displaystyle {\displaystyle \sum }_{i=1}^p}\left(f_i{\rho }_i\right) with 0\le f_i\le 1 \mathrm{and} {\displaystyle {\displaystyle \sum }_{i=1}^p}\left(f_i\right)=1$	[37] (Heinz and Chang, 2001)
Multiple Endmember Spectral Mixture Analysis (MESMA)	$\rho ={\displaystyle {\displaystyle \sum }_{i=1}^p}\left(f_i{\rho }_i\right)+\epsilon with 0\le f_i\le 1 \mathrm{and} {\displaystyle {\displaystyle \sum }_{i=1}^p}\left(f_i\right)=1$ and $RMS={\left(\frac{{{\displaystyle \sum }}_{i=1}^{\lambda }{\left(\epsilon \right)}^2}{p}\right)}^{1/2}$	[39] (Roberts et al., 1998)
Generalized Bilinear Model (GBM)	$\rho ={\displaystyle {\displaystyle \sum }_{i=1}^p}\left(f_i{\rho }_i\right)+{\displaystyle {\displaystyle \sum }_{i=1}^{p-1}}{\displaystyle {\displaystyle \sum }_{j=i+1}^p}\left({\gamma }_{ij}{f}_i{f}_j{\rho }_i\odot {\rho }_j\right)$	[46] (Halimi et al., 2011)
Multi-Linear Model (MLM)	$\rho =\frac{\left(1-P\right){{\displaystyle \sum }}_{i=1}^p{f}_i{\rho }_i}{1-P{{\displaystyle \sum }}_{i=1}^p{f}_i{\rho }_i}$	[34] (Heylen and Scheunders, 2016)

). In MESMA, each endmember was modeled by six spectra (mean spectrum, mean spectrum ± standard deviation, mean spectrum ± 3 × standard deviation and median spectrum). The non-linear contribution was estimated by GBM with its γ_ij parameter and by MLM with the P parameter.

The implementation of the FCLS algorithm comes from the pysptools python’s package. We coded our own version of MESMA based on [39]. GBM and MLM implementation from [44] were used.

3.3. Evaluation Criteria

Three evaluation criteria expressed in percentage of weight abundances (wt %) will assess the performance of the couple spectral preprocessing–unmixing model.

The mean bias (MB) will evaluate the mean error between the estimated abundance ${\widehat{f}}_n$ and the measured one $f_n$ for a pixel n among the region of interest (ROI) image N pixels:

$$MB=\frac{1}{N}\left({{\displaystyle \sum }}_{n=1}^N({\widehat{f}}_n-f_n )\right) ,$$

(3)

The standard deviation around the mean bias (STDB) is as follows:

$$STDB=\sqrt{\frac{1}{N-1}\left({{\displaystyle \sum }}_{n=1}^N{(({\widehat{f}}_n-f_n )-MB)}^2\right)},$$

(4)

The root mean square error (RMSE) is expressed as:

$$RMSE=\sqrt{M{B}^2+STD{B}^2},$$

(5)

4. Results

4.1. Comparison between Reflectance and Preprocessed Spectra

4.1.1. Clay Mixtures

Reflectance spectra for illite–montmorillonite (IM) and montmorillonite–kaolinite (MK) intimate mixtures show that the reflectance level decreased when the montmorillonite abundance increased. For IM mixtures, a spectral shift of 6 nm was observed from I100 to M100 spectra, simultaneously with a slope change in the range 2220–2250nm (Figure 6b). For MK mixtures, the double absorption feature at both 2160 nm and 2200 nm disappeared when montmorillonite abundance was more than 50% (data not shown).

Intrasample spectral variations were more important for reflectance, Hapke and Log(1/R), than for SNV, CR, CWT and 1st SGD whatever IM (Figure 6) and MK mixtures (data not shown). For instance, these variations for reflectance spectra had an amplitude around 2% for illite and 5% for montmorillonite, while for IM mixtures they were around 5% (Figure 6a,b).

The non-linear behavior of intimate mixtures relative to montmorillonite abundance was compared with the linear behavior of aeral ones. For IM mixtures, the mean spectral difference between the two mixture types was not centered to zero, thus highlighting this non-linear behavior (Figure 6). Then, the latter was the most emphasized in clay absorption spectral features for SNV, CR, CWT and 1st SGD suppressing the continuum (Figure 6o,r,u). For reflectance spectra, the non-linearity was quantified up to 1.5% for M20I80 (Figure 6c).

For clay ternary mixtures, same conclusions as before were derived from the analysis of the spectral preprocessing (data not shown).

4.1.2. Montmorillonite–Calcite/Quartz Binary Mixtures

For MC intimate mixtures, the montmorillonite spectral absorption feature started disappearing when calcite abundance was comprised between 80% and 100% (Figure 7b).

The same conclusions as for clay binary mixtures applied for the intrasample spectral variations (Figure 7). For reflectance spectra, they accounted for less than or about 5% (Figure 7b). However contrary to these previous clay mixtures, a better intersample discrimination was noticed with reflectance, Hapke and Log(1/R).

The study of the non-linear behavior of intimate mixtures led to the same conclusions as previously. However, the non-linear effect tended to be higher when calcite abundance dominated, and achieved 4–6% for M20C80 (Figure 7c).

For MQ reflectance spectra, all spectra seemed to overlap for M abundance superior to 11% (Figure 8b). For this reason, the use of spectral preprocessing brought no benefit for spectra separability (Figure 8). MQ aeral mixtures were not shown because of abundances values between aeral and intimate mixtures did not match (cf. Figure 4).

4.2. Performances of Linear Unmixing Methods

4.2.1. Clay Mixtures

Linear unmixing method performances depend on the intrinsic errors due to the unmixing algorithms and the errors due to the non-linearity of intimate mixtures. The first was estimated by applying the unmixing methods on aeral mixtures, and the second would be deduced from the comparison of these results to those obtained with intimate mixtures. These two steps were detailed for clay binary mixtures, while for clay ternary ones; only the global error was given, which accounted for both intrinsic and non-linearity errors (no aeral mixtures available).

For IM aeral mixtures (respectively MK), the best RMSE performance is shown for FCLS applied on SNV, CR, CWT and 1st SGD with RMSE from 2.4% to 3.8% (respectively from 3.2% to 6.2%) (

Table 4. Performances of linear unmixing methods to estimate montmorillonite abundances for illite–montmorillonite (IM) and montmorillonite–kaolinite (MK) aeral and intimate mixtures with different combinations of spectral preprocessings (best results in grey cells with bold font).

Mixture		Aeral								Intimate
		FCLS				MESMA				FCLS				MESMA
		R2	MB	STDB	RMSE	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE
IM	REF	0.66	0.1	15.3	15.3	0.95	−0.2	5.2	5.2	0.14	5.7	22.8	23.5	0.90	7.1	6.9	9.9
	Log(1/R)	0.67	−1.7	15.0	15.1	0.96	−2.3	4.8	5.3	0.23	7.2	21.1	22.3	0.27	9.2	32.4	29.5
	Hapke	0.68	−5.3	14.4	15.3	0.79	−6.5	10.5	12.3	0.30	4.0	19.6	20.0	0.91	0.0	6.6	6.6
	SNV	0.99	2.1	2.2	3.1	0.99	2.2	2.5	3.3	0.98	12.2	3.4	12.7	0.98	12.3	3.4	12.7
	CR	0.98	1.7	3.4	3.8	0.94	1.3	5.3	5.4	0.94	9.4	5.3	10.8	0.95	5.9	5.1	7.8
	CWT	0.99	0.1	2.5	2.5	0.96	0.1	4.6	4.6	0.96	7.1	4.5	8.4	0.80	5.5	10.9	12.2
	1st SGD	0.99	0.1	2.4	2.4	0.97	−0.2	4.1	4.1	0.97	8.3	4.0	9.2	0.97	8.3	4.0	9.3
MK	REF	0.51	−1.3	19.8	19.8	0.97	0.9	3.7	3.8	0.21	12.6	20.9	24.4	0.90	13.3	6.8	15.0
	Log(1/R)	0.53	−2.3	18.9	19.1	0.97	−2.6	3.7	4.5	0.24	15.5	20.2	25.4	0.73	17.4	14.0	21.2
	Hapke	0.58	−3.9	16.9	17.3	0.85	−7.5	9.1	11.8	0.26	11.7	20.1	23.2	0.86	8.8	8.3	12.1
	SNV	0.99	2.1	2.4	3.2	0.99	2.1	2.3	3.1	0.94	14.8	5.4	15.7	0.94	14.8	5.5	15.7
	CR	0.93	1.4	6.0	6.2	0.98	−1.9	3.6	4.1	0.94	14.2	6.4	15.5	0.94	11.3	6.3	13.0
	CWT	0.96	0.0	4.6	4.6	0.99	0.0	2.5	2.5	0.95	12.5	5.2	13.6	0.92	12.7	6.3	14.1
	1st SGD	0.98	0.0	3.5	3.5	0.99	0.3	2.5	2.6	0.94	12.8	5.3	13.9	0.94	12.8	5.3	13.9

). Moreover, CWT and 1st SGD biases were close to zero. Oppositely, FCLS applied on reflectance, Hapke or Log(1/R) led to a large range in M abundance estimates with STDB > 14% (respectively >16%) with an underestimated mean bias (highest bias for Hapke with –5.3%, respectively –3.9%; Figure 9a). For MESMA and both IM and MK aeral mixtures, the best RMSE performance was achieved with whatever reflectance and spectral preprocessing except for Hapke and ranges between 2.6% and 5.3%. However, the range in montmorillonite abundance estimates was strongly reduced compared to FCLS for reflectance, Hapke or Log(1/R) by a maximum factor of 4, but this reduction was not significant for SNV, CR, CWT and 1st SGD (Figure 9b).

For IM intimate mixtures (respectively MK), FCLS methods delivered similar best results with SNV, CR, CWT and 1st SGD with RMSE ranging from 8.4% to 12.7% (respectively from 13.6% to 15.7%; Table 4). Their mean bias in M estimated abundances increase compared to the aeral mixture, from 7.1% to 12.2% for IM intimate mixtures (respectively from 12.5% to 14.8%) and below 2.1% for aeral ones whatever the IM/MK mixtures. For MESMA and IM intimate mixtures, the best results were achieved with reflectance and whatever spectral preprocessing except Log(1/R) with RMSE less than 13% (respectively 15.7% for MK). For IM and MK, a degradation of performances compared to aeral ones was observed around a factor of 2, except for Log(1/R) (degradation by a factor of 5) and Hapke (improvement by a factor of 2 only for IM).

For illite-montmorillonite-kaolinite (IMK) intimate mixtures, better results with FCLS were with SNV, CR, CWT and 1st SGD for RMSE between 9.1% and 14.8% (

Table 5. Performances of linear unmixing methods to estimate montmorillonite abundances for illite-montmorillonite-kaolinite (IMK) intimate mixtures with different combinations of spectral preprocessings (best results in grey cells with bold font).

IMK	FCLS				MESMA
	Montmorillonite
	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE
REF	0.21	0.8	24.3	24.3	0.86	10.8	5.5	12.1
Log(1/R)	0.21	−0.2	23.5	23.5	0.93	8.4	5.1	9.7
Hapke	0.25	−2.9	20.7	20.9	0.77	2.3	6.8	7.1
SNV	0.93	13.9	5.0	14.8	0.92	13.6	5.0	14.5
CR	0.87	8.8	5.2	10.2	0.86	6.3	5.3	8.2
CWT	0.87	7.4	5.2	9.1	0.65	6.0	9.5	11.2
1st SGD	0.89	9.7	4.8	10.8	0.88	10.2	5.1	11.4

). However, the use of SNV led to higher absolute mean biases (MB of around 13.9%) than the others (MB between 7.4% and 9.7%) while the variability in estimates was comparable among them (STDB between 4.8% and 5.2%). For the MESMA method, results performed the same for reflectance and all spectral preprocessing with RMSE ranging from 7.1% (Hapke) to 14.5% (SNV) and similar values of STDB from 5% to 9.5%.

4.2.2. Montmorillonite–Calcite/Quartz Mixtures

For MC aeral mixtures (respectively MQ), linear unmixing provides similar results for reflectance, Log(1/R), CR, CWT and 1st SGD with RMSE from 2.9% to 7% ( respectively from 2.3% to 7.2%;

Table 6. Performances of linear unmixing methods to estimate montmorillonite abundances for MC and MQ aeral and intimate mixtures with different combinations of spectral preprocessings (best results in grey cells with bold font).

Mixture		Aeral								Intimate
		FCLS				MESMA				FCLS				MESMA
		R2	MB	STDB	RMSE	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE
MC	REF	0.96	0.2	4.8	4.8	0.99	0.3	2.5	2.5	0.92	–3.8	8.0	8.8	0.94	7.6	5.4	9.3
	Log(1/R)	0.96	−5.2	4.8	7.0	0.98	−11.5	3.6	12.0	0.93	−8.5	6.7	10.8	0.99	1.6	2.2	2.5
	Hapke	0.91	−20.8	6.6	21.9	0.85	−27.3	8.9	28.7	0.91	−22.6	6.9	23.6	0.87	−23.8	8.3	25.2
	SNV	0.97	12.8	4.1	13.4	0.97	12.9	4.1	13.5	0.86	36.6	12.3	38.6	0.86	36.6	12.3	38.6
	CR	0.98	−5.2	3.6	6.3	0.98	−4.2	4.1	5.9	0.95	19.1	8.6	21.0	0.95	14.3	7.0	16.0
	CWT	0.98	0.0	3.3	3.3	0.98	0.0	3.3	3.3	0.94	21.7	8.3	23.2	0.91	21.4	8.9	23.2
	1st SGD	0.98	0.0	2.9	2.9	0.98	0.1	3.2	3.2	0.96	21.3	7.7	22.6	0.95	22.0	8.1	23.4
MQ	REF	0.93	−0.1	6.2	6.2	0.99	0.4	2.2	2.3	0.29	55.9	17.7	58.6	0.71	53.0	15.1	55.1
	Log(1/R)	0.93	−3.9	6.1	7.2	0.99	−7.2	2.7	7.7	0.29	55.0	17.6	57.8	0.69	53.1	14.7	54.6
	Hapke	0.92	−13.2	6.5	14.7	0.83	−22.4	9.3	24.2	0.32	48.7	18.0	51.9	0.35	48.4	17.8	51.6
	SNV	0.86	35.5	12.3	37.5	0.86	35.5	12.4	37.6	0.00	61.8	20.9	65.2	0.00	61.8	20.9	65.2
	CR	0.98	−5.0	3.4	6.1	0.97	−6.0	3.9	7.2	0.65	51.2	14.4	53.1	0.63	52.1	14.9	54.1
	CWT	0.98	0.0	3.1	3.1	0.98	0.0	3.0	3.0	0.67	53.1	15.6	55.3	0.16	51.2	20.3	55.1
	1st SGD	0.99	0.0	2.8	2.8	0.98	−0.3	3.5	3.5	0.71	53.0	15.3	55.2	0.68	53.1	15.6	55.3

). For MC/MQ mixtures and both FLCS and MESMA, the use of reflectance, CWT and 1st SGD introduced almost no bias in the estimates (≤ 0.4 %).

For MC intimate mixtures, FCLS and MESMA with reflectance or Log(1/R) delivered the best performance with RMSE less than 10.8% and STDB less than 8%. In a least extent, the use of CR with MESMA gave a RMSE of 16% while RMSE accounted for more than 21% for the other spectral preprocessings with both FCLS and MESMA. Furthermore, a large increase in mean bias was observed for FCLS and MESMA coupled to SNV, CR, CWT, Hapke or 1st SGD (MB between 14.3% and 36.6%). Compared to aeral mixtures, poorer performances were obtained with these spectral preprocessing–unmixing method couples by a factor worst than 2, excepted for Log(1/R)–MESMA (no explanation found).

For MQ intimate mixtures, unmixing performance is very low for all linear unmixing methods, with and without spectral preprocessing, with a RMSE higher than 50%. Thus compared to aeral mixtures, performance are the worst for intimate mixtures by a factor of more than 7.

4.3. Performances of Non-Linear Unmixing Methods

4.3.1. Clay Binary Mixtures

For GBM, the best results were obtained with CR, CWT and 1st SGD, leading to RMSE ranging between 13.4% and 13.9% for MK intimate mixtures, and between 7.7% and 9.2% for IM ones (

Table 7. Performances of non-linear unmixing methods to estimate montmorillonite abundances for IM and MK intimate mixtures with different combinations of spectral preprocessings (best results in grey cells bold font).

Mixture		GBM				MLM
		R2	MB	STDB	RMSE	R2	MB	STDB	RMSE
IM	REF	0.56	12.7	14.8	19.6	0.97	9.2	3.8	10.0
	Log(1/R)	0.47	2.0	16.6	16.7	0.97	5.8	4.0	7.0
	SNV	0.98	12.6	3.5	13.0	0.98	11.1	3.3	11.6
	CR	0.95	5.9	4.9	7.7	0.94	11.4	5.9	12.9
	CWT	0.96	7.7	4.4	8.8	0.96	8.2	4.5	9.4
	1st SGD	0.97	8.3	4.0	9.2	0.97	9.4	4.1	10.3
MK	REF	0.56	20.3	14.6	25.0	0.95	13.4	5.1	14.4
	Log(1/R)	0.51	10.0	15.8	18.7	0.95	11.3	5.1	12.4
	SNV	0.94	15.0	5.4	16.0	0.94	13.8	5.3	14.8
	CR	0.94	12.3	6.4	13.9	0.94	12.7	6.1	14.1
	CWT	0.95	12.3	5.2	13.4	0.95	12.4	5.3	13.5
	1st SGD	0.94	12.8	5.3	13.9	0.95	12.9	5.2	14.0

). The worst performance was achieved for reflectance and Log(1/R). Comparing MK with IM mixtures, mean biases were higher with MB of around 12.4% for the first and 7% for the second. The non-linear contribution estimated by GBM was examined through the parameter γ and only for MK mixtures (same conclusions for IM ones; Figure 10). Histograms of γ had narrow peaks close to zero for all spectral preprocessing and reflectance, with CR and 1st SGD having a mean value almost null.

Oppositely, performance of MLM was less dependent on the use or not of spectral preprocessing with RMSE in the range 12.4–14.8% for MK mixtures and 7.0–12.9% for IM ones. Histograms of P, the non-linear contribution of MLM, were centered around zero with a wider distribution covering both positive and negative values (Figure 11). The wider ones were observed for reflectance and Log(1/R).

4.3.2. Montmorillonite–Calcite/Quartz Mixtures

For MC intimate mixtures, GBM combined with reflectance and Log(1/R) gave the best results, respectively with RMSE of 8.2% and 10.8% (

Table 8. Performances of non-linear unmixing methods to estimate montmorillonite abundances for MC and MQ intimate mixtures with different combinations of spectral preprocessings (best results in grey cells bold font).

Mixture		GBM					MLM
		R2	MB	STDB	RMSE		R2	MB	STDB	RMSE
MC	REF	0.92	5.2	6.3	8.2		0.65	26.3	13.4	29.5
	Log(1/R)	0.93	−8.5	6.7	10.8		0.88	12.3	8.1	14.7
	SNV	0.86	35.8	11.9	37.7		0.87	35.3	11.8	37.2
	CR	0.96	16.2	7.8	17.9		0.85	26.3	8.8	27.7
	CWT	0.94	21.6	8.3	23.1		0.88	25.8	8.9	27.3
	1st SGD	0.96	21.3	7.7	22.6		0.89	26.4	8.5	27.7
MQ	REF	0.38	57.4	17.2	59.9		0.70	52.8	15.0	54.9
	Log(1/R)	0.32	54.3	17.2	57.0		0.73	48.8	12.9	50.4
	SNV	0.01	61.8	20.9	65.3		0.01	61.7	20.8	65.1
	CR	0.65	52.2	14.6	54.2		0.00	56.7	23.2	61.2
	CWT	0.67	53.1	15.6	55.4		0.36	53.9	17.4	56.6
	1st SGD	0.71	53.1	15.3	55.2		0.03	51.5	22.5	56.2

). In the least extent, GBM with CR gave a RMSE of 17.9% while using the other spectral preprocessings led to RMSE more than 17.9% due to higher values of MB. Associated histograms of γ were close to zero (Figure 12). MLM results brought more degradation in global performance with RMSE superior to 27%, except for Log(1/R) with RMSE equal to 14.7%. Associated histograms of P were centered to zero for CWT and 1st SGD with narrow peaks, while the histograms of reflectance, CR and Log(1/R) were more spread with negative and positive values (Figure 13). At last for SNV, the histogram was a narrow peak around a mean negative value.

For MQ intimate mixtures, both GBM and MLM presented poor performance with RMSE more than 50% (Table 8).

4.3.3. Clay Ternary Mixtures

For the GBM method, the best results were obtained with CR, CWT and 1st SGD, RMSE respectively of 7.5%, 10.2% and 10.8%, and MB and STDB of the same order of magnitude (

Table 9. Performances of non-linear unmixing methods to estimate montmorillonite abundances for IMK intimate mixtures with different combinations of spectral preprocessings (best results in grey cells with bold font).

IMK.	GBM				MLM
	Montmorillonite
	R2	MB	STDB	RMSE	R2	MB	STDB	RMSE
REF	0.43	12.3	15.3	19.6	0.93	12.0	4.1	12.7
Log(1/R)	0.30	–3.7	18.8	19.2	0.90	8.8	4.6	9.9
Hapke
SNV	0.91	20.5	4.5	21.0	0.93	7.6	5.5	9.3
CR	0.87	5.4	5.2	7.5	0.84	15.6	6.8	17.1
CWT	0.88	8.9	5.1	10.2	0.89	8.3	5.3	9.9
1st SGD	0.89	9.7	4.8	10.8	0.91	10.6	4.9	11.6

). SNV delivered higher absolute mean biases of around 20.5% and RMSE of 21%. For the MLM method, performances were similar whatever the spectral preprocessings (RMSE between 9.9% and 12.7%) except for CR with RMSE of 17.1%. The main errors were driven by MB (15.6%) in comparison with STDB (6.8%).

5. Discussion

5.1. Non-Linearity of Intimate Mixtures Depending on the Mineralogical Composition

Following the increase in montmorillonite abundance, the non-linearity of intimate mixtures was put into evidence compared to aeral mixtures. For two minerals having the same granulometry (around 80 µm), the main contribution of this non-linearity was due to the mean level of reflectance if the minerals have close spectral absorption features such as illite and montmorillonite. Then, with two minerals having no overlapping spectral features such as calcite and montmorillonite, the spectral variations in the absorption features were added as a second contributor to the non-linearity. At last, in the presence of a mineral without any absorption feature such as quartz, the spectral effect of the non-linearity increased drastically. Actually, the reflectance could be decomposed into two components, the volume and the surface reflectances, whose relative contribution depended on the granulometry. First, due to the high reflectance of quartz, the volume reflectance was favored by the multiple scattering inside the mixture and as such increased the chance for a photon to interact with montmorillonite mineral. As a consequence, the non-linearity was observed in the absorption band of montmorillonite, even if its abundance was low as it had already been mentioned by Clark [62]. Second, the surface reflectance depended on the reflective value and the abundance of each mineral. The higher the granulometry, the higher is the surface reflectance relative contribution. Thus for quartz/montmorillonite mixtures, the non-linearity contribution might be partially compensated by the quartz particle size (300 µm).

5.2. Spectral Variability Reduction with Spectral Preprocessings

Application of spectral preprocessings on intimate mixtures spectra showed that 1st SGD, CWT, SNV and CR significantly reduced the intrasample spectral variability while Log(1/R) and Hapke gave an intrasample spectral variability with the same order of magnitude as with reflectance spectra, but sometimes much worse.

Actually, two classes of spectral preprocessings could be distinguished based on their spectral transformation: quasi-linear (1st SGD, CWT) and non-linear (SNV, CR, Log(1/R), Hapke). These two spectral preprocessing classes led to different spectra and introduced some artifacts. CWT decomposes the spectrum into a sum of linear signals and preserves additive properties of the spectra [20] while the 1st SGD is only sensitive to the local spectral variations. SNV is based on a centered-reduced spectrum normalization [17,19] corresponding to a multiplicative correction causing a loss of any additive contribution. Clark et al. [60,62] showed that CR correctly identifies the location of spectral absorption features with spectra having a steep slope continuum but to the detriment of the deformation of the spectral feature [20,62]. Log(1/R) transforms data into the logarithmic space, which, in relationship with Lambert–Beer’s law increases linearity between spectral data and constituent abundance [17,19]. However, additive and multiplicative properties of the spectrum are not suppressed, so that the variability is not reduced. Hapke transforms the data into the albedo space, leading to the same conclusions as Log(1/R).

As a result, the weak intrasample spectral variability of spectra deduced from quasi-linear and non-linear spectral preprocessings 1st SGD, CWT, SNV and CR, was explained by their low sensitivity to the reflectance mean level. On the contrary, the non-linear Hapke and Log(1/R) preserved the intrasample spectral variability of the reflectance. However, the reduction of the latter did not lead to a better discrimination between the different mixture spectra because our samples were homogeneous and thus had low intraclass variability.

5.3. Performance of Linear Unmixing Methods with Spectral Preprocessings

For clay binary aeral mixtures (IM and MK) and reflectance spectra, MESMA is expected to have better performances than FCLS, as it takes into account the intrasample variability of endmember [39]. Improved estimation results were particularly noticeable with reflectance spectra and Log(1/R) having a large intrasample spectral variability, but improvements were very small with the other spectral preprocessings already having a reduced intrasample spectral variability. Overall, FCLS and MESMA had comparable performances with SNV, CR, CWT and 1st SGD. However, SNV and CR carried higher biases in montmorillonite abundance estimation. For calcite and quartz aeral mixtures (MC and MQ), the same conclusions as before could be done in the main lines, with the best results achieved with 1st SGD, CWT and CR.

For clay binary intimate mixtures, the best spectral preprocessings were SNV, CR, CWT and 1st SGD for FCLS whereas for MESMA, Hapke and reflectance spectra were added. MK mixtures had larger montmorillonite estimation errors than for IM ones, oppositely to the clay aeral mixtures. Possible reasons might be that bias values were more important for intimate mixtures than aeral ones. Then for the latter, since there was a strong spectral overlapping of IM mixtures within the range of the intrasample spectral variability (Figure 6a), confusions might occur in finding the best unmixing solution, whereas, confusions might be added due to the non-linear effects for intimate mixtures (more important for MK than IM mixtures). Among the best couples spectral preprocessing/unmixing method, CWT or 1st SGD/FCLS or MESMA were selected. By comparing their results between the aeral and intimate mixtures, the main source of error based on the RMSE came from the non-linearity of the intimate mixtures with an increase by a factor of 2–6. For calcite/montmorillonite intimate mixtures, the best couples spectral preprocessing/unmixing method were REF/FCLS or MESMA. Performances were similar with IM intimate mixtures but better than MK ones. The impact of the non-linearity within the intimate mixtures compared to the aeral mixtures increased by a factor of 2–4. For quartz/montmorillonite intimate mixtures, the non-linearity effect due to the presence of quartz was too strong to be taken into account by the linear unmixing methods and the performances were very low. For ternary clay intimate mixtures, CR or CWT/FCLS or /MESMA were the best couples and their performance had the same order of magnitude as those of clay binary mixtures.

As a first conclusion, when two minerals present overlapped absorption features (clay mixtures), the best spectral preprocessings were those enhancing the slight spectral local changes, such as SNV, CR, CWT and 1st SGD. However, SNV was not recommended since it might bring higher errors on mineral abundance estimation in some cases. Oppositely, when no overlap occurred (mixtures with calcite), Hapke, and reflectance spectra performed the best because they rely on the global variation of the spectrum (i.e., continuum), which varied linearly with the addition of minerals.

5.4. Performance of Allunmixing Methods with Spectral Preprocessings

For clay binary intimate mixtures and whatever linear or non-linear unmixing method, results of montmorillonite abundance estimation were in the range 10.8–25.4% (RMSE) for MK mixtures and 7.0–29.5% for IM ones. MESMA, FCLS and GBM had comparable performances when using SNV, CR, CWT and 1st SGD, the same as MLM combined with whatever spectral preprocessing. MLM algorithm might compensate the higher intrasample spectral variability with the P non-linear parameter, because this spectral variability had a larger amplitude than the non-linearity of the intimate mixture for similar grain size in clay mixtures [34]. Revel et al. [63] have also observed the prevalence of intraclass variability over non-linear effects on the performance of unmixing methods in another context.

For calcite/montmorillonite intimate mixtures, RMSE ranged between 8.2% and 38.6%. GBM had a very close performance with reflectance and Log(1/R) compared to FLCS and MESMA with reflectance, while MLM performed the worst. Compared to clay binary intimate mixtures, this last case may be explained because the P parameter of MLM had negative values. Heylen et al. [34] mentioned that the use of MLM could lead to errors in the case of high reflective materials, like calcite.

For quartz/montmorillonite intimate mixtures, results of montmorillonite abundance estimation were very poor with RMSE around 50.4–65.3%. The difficulties raised in this study with mixtures containing quartz are also noticed by Asadzadeh et al., Viscarra Rossel et al. and Mulder et al. [15,26,29] since detecting and quantifying quartz are still major limitations for soil mineralogy spectroscopy in the range 400–2500 nm. Indeed, some studies chose to neglect the impact of quartz abundance [26]. However, Debba et al. [16] demonstrated that taking into account quartz abundance as an endmember for unmixing produces more accurate results.

For clay ternary intimate mixtures, besides the increasing impact of non-linear effects from clay binary to clay ternary intimate mixtures, the performance was globally similar whatever the unmixing method with montmorillonite abundance estimation in the range 7.5–24.3% (RMSE). The best performance was obtained for CWT or 1st SGD/FCLS or MESMA or GBM or MLM.

Globally several trends could be deduced, except for the quartz/montmorillonite mixtures:

• For clay intimate mixtures either binary or ternary and whatever unmixing methods, the best spectral preprocessings were CWT and 1st SGD while for calcite/montmorillonite intimate mixtures, the best results were obtained without the use of spectral preprocessing. Consequently, it was recommended to use quasi-linear spectral preprocessings (CWT or 1st SGD) when the absorption bands of minerals overlap and keep the reflectance in the other cases (calcite). Attention should be paid with the use of SNV leading to higher biases in montmorillonite abundance estimations, and Log(1/R)/Hapke leading to non robust results dependent on the mixture type and the unmixing method. The only best results for Hapke were obtained with MESMA for clay intimate mixtures (IM: RMSE of 6.6%, MK: RMSE of 12.1% and IMK: RMSE of 7.1%) with minerals of similar granulometry (80 µm), which is in agreement with Heylen and Scheunders [31] (RMSE less than 2% for alunite–quartz mixtures having similar grain size). However, with minerals of different granulometry (calcite: 70 µm, quartz: 300 µm), performance was poor because it violated the main assumption of the Hapke algorithm used in this paper [44].
• Similar performances in terms of RMSE were noticed between the use of linear and non-linear unmixing methods. As a result, the use of the simplest linear unmixing method, FCLS, was advised, coupled with the simplest spectral preprocessing, 1^st SGD, for clay intimate mixtures and without the use of spectral preprocessing with calcite/montmorillonite.
• The error in montmorillonite abundance estimation achieved for the best couple mentioned in Section 2, a RMSE of 9.2% for IM mixtures, 13.9% in MK mixtures, 10.8% in clay ternary mixtures and 8.8% for MC mixtures. These results were better than those obtained by [12] using a geometrical analysis for montmorillonite–illite–kaolinite mixtures (RMSE 15.5%). For more complex mixtures, the performance gave an RMSE of around 8% using a regression tree for smectite–kaolinite–muscovite–calcite–quartz mixtures [26] and an RMSE of around 3.4% using a multivariate analysis for smectite–illlite–kaolinite–carbonate–quartz mixtures [29].

6. Conclusions

A comparative study was carried out in order to assess the performance of combining unmixing methods (two linear and two non-linear methods, FCLS and MESMA, GBM and MLM, respectively), with and without the use of six spectral preprocessings (SNV, CR, CWT, Hapke, 1st SGD and Log(1/R)), for the estimation of montmorillonite abundance. The objective of this work was to analyze the sources of non-linearities as a function of the mixture type (aeral versus intimate), the number of minerals (binary and ternary) and their nature (clay, calcite and quartz).

The major results included the following: (i) spectral preprocessing SNV, CR, CWT and 1st SGD reduced the spectral intrasample variability, (ii) the benefit of the spectral preprocessings CWT and 1st SGD occurred when spectral absorption features of minerals overlapped, whereas the reflectance spectra without spectral preprocessing performed the best when no overlap occurred, (iii) SNV carried non-linear effects, which led to biases for montmorillonite estimation and the use of Log(1/R) and Hapke sometimes led to non robust results, (iv) unmixing on mixtures with quartz achieved the worst performance with RMSE higher than 50%, (v) linear and non-linear unmixing methods had similar performance and so the use of FCLS, the simplest method, might be recommended and (vi) the most robust couple of spectral preprocessing and unmixing method was 1st SGD and FCLS for the clay binary mixtures and reflectance and FCLS for the mixtures with calcite, with RMSE ranging between 8.8% and 13.9%. For soil swell–shrinking application, our results with ternary and binary clay and clay–calcite mixtures gave an RMSE less than 15%, the use of unmixing methods could be of interest in order to improve the classification of expansive soils.

This study was only carried out with mineral mixtures manually generated in laboratory conditions, in order to study non-linear interactions due to mineralogy alone. However, these mixtures did not represent real soil composition, mainly because they did not include the presence of soil organic matter and water content and because of the homogeneous mineral content (little variability in a mixture, with a controlled grain size). Thus, future work will be to carry out the same comparative study applied on real soils with very high resolution hyperspectral imaging in order to take into account the impact of soil composition like mineralogy and, water content and organic matter content [7,29], as well as environmental factors like soil surface roughness [64] and shadow effects [64,65,66,67] and the presence of either sparse green or dry vegetation with soil within a pixel [49,68,69,70]. As quartz and carbonates minerals possess distinct features in long wave infrared (LWIR, 7.5–14 µm) [56,71], we propose in a further perspective to use the combination of SWIR and LWIR spectral domains in order to improve clay estimation in soil.