Evaluating Endmember and Band Selection Techniques for Multiple Endmember Spectral Mixture Analysis using Post-Fire Imaging Spectroscopy

Abstract

Fire impacts many vegetated ecosystems across the world. The severity of a fire is major component in determining post-fire effects, including soil erosion, trace gas emissions, and the trajectory of recovery. In this study, we used imaging spectroscopy data combined with Multiple Endmember Spectral Mixture Analysis (MESMA), a form of spectral mixture analysis that accounts for endmember variability, to map fire severity of the 2013 Rim Fire. We evaluated four endmember selection approaches: Iterative Endmember Selection (IES), count-based within endmember class (In-CoB), Endmember Average Root Mean Squared Error (EAR), and Minimum Average Spectral Angle (MASA). To reduce the dimensionality of the imaging spectroscopy data we used uncorrelated Stable Zone Unmixing (uSZU). Fractional cover maps derived from MESMA were validated using two approaches: (1) manual interpretation of fine spatial resolution WorldView-2 imagery; and (2) ground plots measuring the Geo Composite Burn Index (GeoCBI) and the percentage of co-dominant and dominant trees with green, brown, and black needles. Comparison to reference data demonstrated fairly high correlation for green vegetation and char fractions (r² values as high as 0.741 for the MESMA ash fractions compared to classified WorldView-2 imagery and as high as 0.841 for green vegetation fractions). The combination of uSZU band selection and In-CoB endmember selection had the best trade-off between accuracy and computational efficiency. This study demonstrated that detailed fire severity retrievals based on imaging spectroscopy can be optimized using techniques that would be viable also in a satellite-based imaging spectrometer.

1. Introduction

Fire behavior, size, and severity are changing in the western United States [1,2,3]. To fully comprehend these changes, techniques to reliably map fire effects over large areas are required. The most common means of assessing fire severity with remote sensing is using the Normalized Burn Ratio (NBR) and its derivatives [4,5]. NBR based techniques are popular because of their simplicity, ease of use, and relatively good performance when compared to field measurements [6,7]. The NBR, however, is sensitive to variations in soil brightness [8], vegetation type [9], and vegetation density [10]. These limitations constrain the use of NBR-based indices for consistent fire severity assessments over large areas within heterogeneous ecosystems [10]. Variation in fire severity can effectively be broken down into detectable differences in the relative abundance of char, green vegetation, dead vegetation, and bare soil; remotely sensed fire severity assessments are therefore essentially based on mixtures composed of these four constituents. Under this paradigm, Spectral Mixture Analysis (SMA), in which reflectance is assumed to a be a linear combination of components or endmembers at a subpixel level [11,12], represents a viable alternative to NBR-based analysis, potentially overcoming NBR’s sensitivity to different cover type and soil brightness variation [13,14]. SMA has been used previously to characterize tree mortality and soil char cover [10,15,16,17,18].

Lentile et al. [19] defines fire severity as the immediate ecosystem impacts of fire; burn severity, in contrast, is defined as the combined shorter and longer term ecosystem impacts and response to fire. Key and Benson [4] also separates first-order effects, which are related to the fire only, and second-order effects, which are related to other environmental (e.g., wind and rain) and vegetative processes. Veraverbeke et al. [20] followed these conventions and used the term fire severity when images were acquired soon after the fire (e.g., within the first month) and second-order effects can safely be neglected.

In SMA, an endmember class is one of the categories the image will be divided into, and an endmember is the individual training spectra that are constituents of each endmember class. In simple SMA, only one endmember represents each endmember class across the image. Multiple Endmember SMA (MESMA), in contrast, allows the endmember representing each endmember class to vary on a per-pixel basis [21]. This approach accounts for the variability that may exist within each endmember class and further allows for consistency and accuracy across ecosystems where there could be considerable spectral variability within an endmember class [21,22].

SMA is particularly well suited for imaging spectrometry data, as the large number of bands provide additional information for cover type discrimination in critical, but spectrally narrow regions; for example, the red edge for green vegetation discrimination [23]. The Airborne Visible/Infrared Imaging Spectrometer (AVIRIS), samples spectra from 350 nm to 2500 nm at a nominal sampling interval of 10 nm [24]. There is currently no spaceborne imaging spectrometer with a similar signal-to-noise ratio to AVIRIS and global coverage with frequent return intervals, which currently impedes the use of imaging spectroscopy for fire severity comparisons at regional scales and over longer time periods. However, there are several proposed spaceborne imaging spectrometers, including the Environmental Mapping and Analysis Program (EnMAP, [25]), the PRecursore IperSpettrale della Missione Applicativa (PRISMA, [26]), and the Hyperspectral Infrared Imager (HyspIRI, [27]), that would make using imaging spectroscopy to monitor fire severity possible at regional to global scales within the decade.

A drawback of MESMA, particularly if three or more endmember classes are allowed to be modeled within a single pixel, is that the number of endmember combinations tested can be large. A variety of approaches have been developed that attempt to select the minimum number of spectra in a library needed to represent within endmember class variation and therefore eliminate redundant spectra. Many techniques focus on automated or semi-automated endmember selection, which identifies pure pixels through the extremes of image data [28] or the construction of synthetic endmembers based on image data [29]. Other techniques first create a large spectral library from various sources and then reduce its size to achieve a library that is both parsimonious and captures each endmember class’s variability [30,31,32]. Several means exist to evaluate which endmembers to keep. Generally, criteria evaluated focus on either which endmembers best represent their endmember class [33,34,35], or which endmembers best model the library as a whole [31]. All techniques have slightly difference balances between capturing spectral variability and creating efficient libraries. While there are several studies reviewing endmember extraction techniques [36,37,38], the comparison of endmember selection techniques for MESMA is rare [32], and to our knowledge, there has been no evaluation of various endmember selections techniques for cover fraction identification (MESMA with multiple non-shade endmembers allowed per pixel).

Individual bands in imaging spectroscopy tend to be highly correlated, and the inclusion of all bands in image analysis techniques increases computational times and can decrease accuracy [36,39]. Data reduction techniques are therefore often applied imaging spectroscopy data sets prior to SMA [36]; two of the more common are Principal Component Analysis (PCA; [40]) and maximum noise fraction (MNF; [41]). PCA, MNF, and similar techniques reduce data dimensionality based on the spectral properties of the image; however, they do not necessarily do so in a manner that improves separability between endmember classes [42]. Asner and Lobell [43] proposed a data reduction technique designed specifically to improve the accuracy of SMA when applied to plant cover with AutoSWIR. AutoSWIR uses a priori optical properties of leaf, litter, and soil material to select critical spectral regions to use for spectral unmixing. Somers et al. [44] proposed a technique similar to AutoSWIR in stable zone unmixing (SZU), which uses variability within and between endmember classes to select spectral regions to use in SMA. Since SZU is based solely on the input spectral library, it has the built-in assumption that the endmembers in a spectral library are representative of the variability of the spectra in the image to be analyzed; however, it has the advantage of selecting spectral regions specific to the problem. SZU has been demonstrated to be effective in invasive species monitoring, soil type classification, and oil spill detection and tracking [44,45,46]. Neither autoSWIR nor SZU specifically addresses the highly correlated nature of adjacent bands. Somers and Asner [47] proposed a further refinement on SZU, uncorrelated SZU (uSZU), which selects bands that capture the maximum variability within and between endmember classes while eliminating highly correlated bands. uSZU was shown to have improved cover abundance estimate accuracy and performance times compared to SZU [47].

In this study, we test the ability of imaging spectroscopy based MESMA to derive cover fractions that correspond with field measurements. Additionally, we test the specific effects of four different endmember selection techniques and one band reduction technique on MESMA’s accuracy and processing time in the context of deriving indicators of fire severity of a large wildfire using imaging spectroscopy.

2. Methodology

In order to test the effectiveness of derived cover fractions at modeling the actual mixed composition of pixels, four post-fire AVIRIS flight lines were processed using MESMA. We used four endmember selection techniques and tested each technique with both the full AVIRIS spectra and a reduced spectral subset determined from a band selection algorithm. The result was eight separate spectral libraries. Each of these eight libraries was used to perform a separate run of the MESMA algorithm (Figure 1). In order to evaluate MESMA derived fractions correspondence with fire severity, two separate validation data sets were used. First, the various cover fractions generated by the eight MESMA runs were compared through linear regression against randomly selected and manually classified plots on WorldView-2 scenes. Using MESMA derived fractions from a spectral library that performed well, the relationship of the cover fraction with two field measurements, Geo Composite Burn Index (GeoCBI), and the percent cover of green, brown, and black trees over an area, was evaluated.

2.1. Study Area

The Rim Fire occurred in California’s Sierra Nevada, starting in the Stanislaus National Forest and burning into Yosemite National Park (Figure 2). The Rim Fire started on 17 August, 2013 and was contained on 24 October 2013. According to the official fire perimeter statistics from the state of California (http://frap.fire.ca.gov/) the Rim Fire burned 104,131 ha and was the largest fire ever recorded in the Sierra Nevada and the fourth largest in California [3]. Several vegetation types were affected by the fire: coniferous forests dominated by Ponderosa pine (Pinus ponderosa), Lodgepole pine (Pinus contorta), red fir (Abies magnifica), and mixed conifer stands; deciduous and evergreen hardwood forests made up of various lower montane species including blue oak (Quercus douglasii), valley oak (Quercus lobata), interior live oak (Quercus wislizeni), and California buckeye (Aesculus Californica); annual grass; and shrubs, primarily composed of many common upper and lower montane Sierra Nevada chaparral species such as manzanita (Arctostaphylos spp.) and California lilac (Ceanothus spp.).

2.2. AVIRIS Imagery and MESMA

2.2.1. Preprocessing

Level 2 ortho-reflectance AVIRIS flight lines were acquired from the Jet Propulsion Lab (JPL) ftp site (https://aviris.jpl.nasa.gov/data/AV_HyspIRI_Prep_Data.html). The flight lines used in this analysis were f131117t01p00r07, f131117t01p00r08, f131117t01p00r09, and f131117t01p00r10, all acquired on 17 November 2013. Before retrieval, the images had been orthorectified. Surface reflectance was retrieved using ATmospheric REMoval program (ATREM; [48]). Images were further manually georeferenced using United States Geological Survey (USGS) 2004 1m digital orthophoto quarter quad (DOQQ) maps (http://nationalmap.gov/viewer.html) as base maps. Images were warped using a first-degree polynomial with nearest neighbor resampling. Images were resampled to a uniform 14.8 m on a side pixel size. The flight lines cover approximately 93% of the Rim Fire burn scar. Wavelengths 365–405 nm, 1325–1420 nm, 1810–2040 nm, and 2450–2500 nm were known to either be associated with atmospheric water absorption features or found to contain significant atmospheric artifacts and a low signal-to-noise ratio from visual inspection; bands within these wavelengths were not considered in further processing, leaving 177 AVIRIS bands for use in the analysis. In addition, flight lines over the same area from 26 June 2013 were used as sources for the spectral library used for the MESMA process. The flight lines used for this purpose were f130626t01p00r13, f130626t01p00r14, and f130626t01p00r15. These flight lines were processed using procedures already described.

2.2.2. Spectral Library

Five different endmember classes were used to unmix the Rim Fire images: green vegetation (GV); non-photosynthetic vegetation (NPV), which includes dead needles, forest litter, bark senesced grass and other non-photosynthetically active plant material; soil, which includes bare soil and rock; char, which includes charred plant material and ash; and shade. Typical sources for endmembers include either spectra derived from the field, laboratory, or spectra derived from the imagery itself using pure pixels [49].

In this study, a combination of sources were used (

Table 1. The number and origin of endmember used for the source spectral library. Analytic Spectral Devices is abbreviated as (ASD) and denotes field collected samples spectrally analyzed at a lab. Images acquired through Airborne Visible/Infrared Imaging Spectrometer are abbreviated as (AVIRIS). Spectra developed during a study in Wind River Research Forest [51] are denoted as Wind River. Green vegetation is abbreviated as (GV) and non-photosynthetic vegetation as (NPV).

	Char	GV	NPV	Soil	Total
17 November AVIRIS	457	308	0	358	1227
26 July AVIRIS	0	1739	245	510	2634
ASD	21	0	3	46	70
Wind River	0	498	129	139	766
Total	478	2545	377	1053	4453

). Field samples of charcoal, non-photosynthetic vegetation and substrate were collected and analyzed at the JPL spectroscopy lab [50]. Measurements in the 400–2500 nm spectral domain were obtained with an Analytic Spectral Devices (ASD) Full Range Spectrometer under artificial lighting conditions with an ASD Pro lamp. Reflectance was calibrated using a white Spectralon panel. We also used spectra collected from AVIRIS imagery, acquired on 26 June 2013 and 17 November 2013. To extract these spectra, polygons of areas composed solely of one of the four endmember classes (not including shade) were identified in the images. The spectra of all pixels within each polygon were extracted and labeled by endmember class. The spectral profile of each pixel was then manually reviewed to assure that it was consistent with the endmember class it was labeled as. To capture spatial scales that were not otherwise available in our analysis, an existing spectral library created from Wind River Research Forest (in Southern Washington state) was used [51]. This library included collections of leaves (branches) and tree bark, measured with an ASD in the field. Tree crowns measured in situ via a canopy crane were also included in this library. Although the spectra in the Wind River library were acquired in a different geographic location, many of the species and genera are common for both locations. Although Roberts et al. [51] also included stand level spectra measured with AVIRIS, these spectra were not included in our study.

2.2.3. Band Selection

Once the spectral library was finalized, the uSZU band selection technique [47] was performed to identify bands that optimize endmember class discrimination. uSZU is an automated band selection technique which attempts to select bands which maximize variation between endmember classes in a spectral library while minimizing the correlation between selected bands. uSZU assigns each band a Stability Index (SI) value based on the interclass variability divided by the intraclass variability for each class [47]:

$$S{I}_i=\frac{\Delta inter,i}{\Delta intra,i}=\frac{\left|{{\displaystyle \sum }}_{j=1}^k{R}_{mean,j,i}\right|}{1.96\times {{\displaystyle \sum }}_{j=1}^k{\sigma }_{j,i}}$$

(1)

where k is the total number endmembers, R_mean,j,i is the mean reflectance for endmember class j at wavelength i, and σ_j,i is the standard deviation of class j for wavelength i. The band with the highest SI value is selected. Then a spectral correlation value (Corr) is calculated:

$$Corr\left(X,Y\right)=\frac{cov\left(X,Y\right)}{{\sigma }_x{\sigma }_y}$$

(2)

where cov(X,Y) is the covariance between the selected waveband, X and each remaining waveband, Y; and σ is the standard deviation. All bands with corr values above a predetermined threshold c are then discarded. The process is repeated, with the threshold for correlation needed to discard a band decreasing by the value i with each iteration. For this study, values of c = 0.99 and I = 0.001 were used, both of these values were tested in Somers and Asner [47] and found to give acceptable accuracy. MESMA that was based on uSZU band reduced spectral libraries will be noted by including the term “uSZU” in the name.

2.2.4. Endmember Selection

Several approaches have been developed for determining the relative value of individual endmembers for representing their endmember class. One method for doing this is a count-based (CoB) approach, where each endmember is selected iteratively, and, using simple SMA, tested to see how many endmembers it can successfully model within the spectral library, using a predetermined Root Mean Square Error (RMSE) threshold to define success [52]. When this approach is applied to an endmember within its own endmember class, it is called (In-CoB); a desirable endmember will have a high In-CoB number (indicating a large number of other endmembers within the endmember class are derivative). This technique can also be applied against all the endmember classes which the endmember does not belong (Out-CoB), in this case, a desirable endmember will have a low Out-CoB number, indicating that this endmember would minimize confusion with other endmember classes.

Another approach for selecting endmembers is Endmember Average Root Mean Square Error (EAR), which evaluates each endmember’s ability to model all other within class endmembers based on a summed RMSE [33]. EAR is calculated using the following formula:

$$EA{R}_{A_i}=\frac{{{\displaystyle \sum }}_{j=1}^{n}RMS{E}_{A_i{A}_j}}{n-1}$$

(3)

where A is an endmember class, A_i is the selected single endmember, and A_j are each of the other endmembers within the endmember class, and n is the total number of spectra in class A. A smaller EAR value is more desirable. Another approach, Minimum Average Spectral Angle (MASA, [34]), is similar to EAR but evaluates the summed spectral angle [53] instead of the RMSE fit. Both EAR and MASA evaluate only within an endmember class and do not evaluate interclass confusion.

Two techniques based on EAR, MASA, and CoB were used in this study. The first technique used a combination of EAR, MASA, and CoB to select three endmembers, it will therefore be abbreviated to “EMC”. The combined EMC technique selected three endmembers for each endmember class: an endmember which minimizes the EAR value, an endmember which minimizes the MASA value, and an endmember which maximizes the In-CoB value. If multiple endmembers had the same In-CoB value, the endmember with the smallest Out-CoB value was selected. If the same endmember is selected through multiple EMC criteria (for example, the spectra that minimizes MASA also minimizes EAR), then fewer than three endmembers were used for that endmember class. The In-CoB [35] technique was also used in this study. For this technique, the In-CoB value of each endmember in an endmember class was evaluated, and any endmember with a unique In-CoB value was selected. If multiple endmembers have the same In-CoB value, the endmember with the minimum EAR value was selected. Another approach to endmember selection is Iterative Endmember Selection (IES) [31], which first picks two endmembers that maximize the performance of two-endmember SMA classifying the entire spectral library, as determined by using Cohen’s kappa [54]. The algorithm then iteratively adds and removes endmembers in order to maximize kappa until the smallest spectral library that maximizes kappa is developed. While the other endmember selection approaches focus on within-class variability, IES explicitly considers confusion of endmembers between classes. In this study, the parameter of 0.025 was used as an RMSE constraint for the two-endmember SMA needed to classify the spectral library [32]. Spectral libraries processed in this way will be referred to as “IES” for the remainder of this study.

IES generally results in relatively large spectral libraries that can make the resultant MESMA analysis computationally expensive. Roberts et al. [55] proposed a method for reducing the size of IES generated spectral libraries through an iterative process, wherein the brightest endmember in each endmember class is initially chosen and all endmembers are modeled as a mixture of that endmember and shade are subsequently removed from the spectral library. Following the first pass, a new set of bright endmembers are selected from the reduced library, targeting endmembers that are most spectrally distinct from the first set. The process is repeated iteratively until either there are no more endmembers within a class in the reduced library, or new endmember selections fail to further reduce the library. The goal of this process is to identify the smallest set of endmembers for each class that also fully characterizes the spectral diversity of that class. The final result is a spectral library that retains the high accuracy of IES, while significantly reducing the number of endmember combinations. Spectral libraries processed in this way will be referred to as “Reduced IES”.

2.2.5. MESMA

Four endmember MESMA (a maximum mixture of three endmember classes and shade could be modeled for any given pixel) was performed on the post-fire images using the resultant eight spectral libraries. In order for a model with more endmembers to be chosen over a model with fewer endmembers, the model with more endmembers must improve the model’s RMSE for the pixel by a minimum threshold of 0.007 [35]. The endmember and shade fraction were constrained to be physically reasonable (no fractions below 0 or above 1 were allowed). The shade fraction was further constrained to be no higher than 0.8. A maximum allowed RMSE of 0.025 was also used to constrain the MESMA model [21].

2.3. Validation

2.3.1. WorldView-2 Imagery

High spatial resolution imagery provides a validation source for cover fractions over large portions of the fire area. WorldView-2 is a commercially operated very high spatial resolution, multispectral satellite imager. Twelve WorldView-2 images from 28 September 2013 that together captured about 97% of the Rim Fire burn scar (Figure 2) were purchased from Satellite Imaging Corp. The images had a pixel resolution of 2 m. After acquisition, images were converted into reflectance using Atmospheric and Topographic CORrection (ATCOR) as implemented in the PCI Geomatica software package (PCI Geomatics, Geomatica version 2014) using a LiDAR data set acquired by the National Center for Airborne Laser Mapping (NCALM-University of Houston).

In order to compare the AVIRIS-derived MESMA analysis and the WorldView-2 images, 120 square samples were randomly selected. The samples were 119.2 m on a side, or equivalent to an eight by eight square of AVIRIS pixels, this was chosen to reduce problems of differences in georeferencing between the images, as well as reduce the effect of pixels which are only partially within the sample area in the AVIRIS imagery. Areas within the sample squares were manually interpreted and digitized from the WorldView-2 imagery into five cover classes: GV, NPV, char, soil, and shadow (see Figure 3). The spectral resolution of WorldView-2 typically allowed for easy discrimination of GV. The NPV, char, and soil classes were not always easy separable spectrally, but contextual clues such as shape and texture were often helpful in separating these classes. There were likely some pixels that were misclassified, but overall this was considered a minor source of error.

For both the eight MESMA unmixing results and the manually interpreted high-resolution imagery, the average cover for each cover class was calculated for each point, then the fractions were shade normalized, in which the shade cover fraction was removed by dividing each endmember by the total percent cover of all non-shade endmembers [56]. This allowed the relative abundance of the non-shade endmembers classes to be compared by partly correcting for the effects of differential solar geometry between the two images [13]. The quality of a linear fit for the relationship between each of the eight MESMA analysis and the manually identified imagery based on validation plots was assessed.

2.3.2. Field Plots

Thirty-four field plots of fire severity were collected in October 2013. The plots were located in the southeastern part of the fire and were, due to permitting restrictions, all taken within Yosemite National Park. Five severity was assessed in 30 by 30 m square shaped plots. The plot’s center coordinates were recorded with a handheld Trimble GeoExplorer 6000 series GeoXH GPS device (these unit average approximately 1-m x and y error, see Figure 2 for plot locations). The plots were located no less than 200 m apart and were taken in areas of relatively homogeneous fuel type and fire severity. Four of the 34 original plots were removed from the analysis since they were located outside of AVIRIS images of 17 November 2013.

GeoCBI [57], a modified version of Composite Burn Index (CBI) [4], was measured in the field plots. GeoCBI divides the ecosystem into five different strata: (i) substrates; (ii) herbs, low shrubs and trees less than 1 m; (iii) tall shrubs and trees of 1 to 5 m; (iv) intermediate trees of 5 to 20 m; and (v) tall trees taller than 20 m. In the field form, 20 different factors can be rated (e.g., soil and rock cover/color change, percent leaf area index change, char height) but only those factors present and able to be reliably evaluated for that strata are considered. Ratings are given on a continuous scale between zero and three and the results are averaged per stratum. These stratum averages are weighted according to their fractional cover within each plot to compute a final GeoCBI rating. GeoCBI scores range from zero to three, with three representing the highest severity. In addition to GeoCBI, every dominant and co-dominant tree within each plot was identified. Trees were then assigned to three possible classes based on their needle color: green, brown, or black.

To compare the 30 by 30 m GeoCBI plots with the 14.8 m pixel size, a weighted plot level averaging of pixels was performed. Weights were assigned based on each pixel’s percentage of area within a plot. A weighted average of all pixels within an individual plot was then calculated, and compared with the plot level data.

3. Results

3.1. uSZU Band Selection

After uSZU of the 177 original bands, 20 were retained in the analysis. Selected bands ranged across the spectrum with a preference for the red edge (five bands were selected between 670 and 735 nm). uSZU also appears to have selected several bands that are at the edges of wavelengths removed due to atmospheric noise, suggesting the sensitivity of the technique to spectral artifacts. Figure 4 shows the final selection of bands using the uSZU algorithm.

3.2. Endmember Selection and Processing Times

The number of endmembers selected by each technique varied considerably (

Table 2. The number of endmembers and modeling times of various combinations of endmember and band selection techniques. Modeling time was assessed on a 640 × 530 pixel Airborne Visible/Infrared Imaging Spectrometer (AVIRIS) scene within the Rim Fire. Models were run using a computer with an Intel Xeon central processing unit E7-4850v2 at 2.3 GHz with 48 cores and 512 gigabytes of random-access memory. Models were run in a manner to minimize run time with computer resources. Uncorrelated Stable Zone Unmixing (uSZU) represents a feature reduction technique. EMC, In-COB, and IES are all endmember reduction technique for Multiple Endmember Spectral Mixture Analysis (MESMA). Total models represents the number of models need to run four endmember MESMA.

Class	EMC	uSZU EMC	In-CoB	uSZU In-CoB	IES	uSZU IES	Reduced IES	uSZU Reduced IES
Char	2	3	5	5	10	10	5	7
GV	2	2	14	17	31	42	25	40
NPV	3	2	11	11	36	40	32	37
Soil	2	2	15	16	55	53	35	38
Total Models	83	83	5729	7071	115,543	156,821	45,327	92,415
Processing Time	0.98	0.44	14.88	8.70	151.85	66.76	78.36	54.45

). EMC has an inherent limit to the number of endmembers that may be selected; in EMC no more than the number of endmember classes multiplied by three endmembers are kept (the maximum possible number of endmembers in an EMC reduced library is therefore 12 in this study). uSZU generally had the effect of increasing the number of endmembers selected for all techniques except EMC. The increase was generally moderate, resulting in a 10–25% increase in the number of endmembers. Despite the larger number of endmembers selected, the processing time for the MESMA algorithm was dramatically reduced for all band reduction techniques, with time reductions ranging between 30% and 55%.

3.3. Unmixed Images and Overall Model Comparison

The eight endmember libraries produced broadly similar patterns (Figure 5). The EMC method without band selection stands out as modeling more NPV than the others, notably modeling NPV in areas that other approaches modeled as GV. Unmixed images generated from libraries without uSZU band selection appear to have slightly more NPV than their uSZU counterparts, but otherwise give similar results.

While the models were qualitatively similar, there were quantifiable differences between them (

Table 3. Percent of total pixels within the Rim Fire boundary that were successfully modeled and the percent of pixels containing each endmember class for various endmember selection techniques. In total, approximately 7,850,000 pixels within the Rim Fire boundary were considered. A single pixel was allowed to contain up to three endmember classes (not including shade).

	EMC	uSZU EMC	In-CoB	uSZU In-CoB	IES	uSZU IES	Reduced IES	uSZU Reduced IES
Modeled	87.0%	81.7%	83.2%	79.1%	93.0%	92.0%	85.7%	86.4%
Char	36.9%	53.1%	47.3%	36.9%	33.4%	42.9%	41.6%	34.4%
GV	54.1%	52.9%	43.1%	39.9%	35.3%	39.4%	36.8%	38.1%
NPV	38.8%	21.9%	27.7%	22.0%	23.5%	21.2%	30.0%	21.65
Soil	7.9%	7.5%	11.5%	26.2%	41.8%	28.6%	17.5%	28.4%

). The number of pixels successfully modeled, defined by the number of pixels where any endmember combination tested by the MESMA algorithm resulted in an RMSE below the threshold of 0.025, ranged from 79.1% of the fire area (uSZU In-CoB) to 93.0% of the fire area (IES). Soil was the most variable endmember class between the models, with one scene containing as few as 7.5% of pixels with any modeled soil cover (uSZU EMC) to 41.8% of pixels (IES). Char and GV were the most commonly modeled endmember classes with between 36.9% and 53.1% of pixels containing char and between 35.3% and 54.1% of pixels containing GV.

3.4. Endmember Sources in Model Selection and the Image

Most spectral libraries selected endmembers from all sources (

Table 4. The source of endmembers for each endmember reduction technique.

Spectra Source	EMC	uSZU EMC	In-CoB	uSZU In-CoB	IES	uSZU IES	Reduced IES	uSZU Reduced IES
AVIRIS	9	9	15	18	41	45	20	30
ASD	0	0	5	4	12	14	9	12
Wind River	0	0	25	27	79	86	68	80

). In the initial spectral library, a majority of spectra were derived from AVIRIS imagery, however, a disproportionally small amount of spectra from this source were selected by most endmember reduction techniques. This suggests that, by the metrics used by the endmember reduction techniques, most image derived spectra were redundant and were therefore eliminated. In contrast, spectra created in Wind River [51] were disproportionally retained by most endmember selection techniques. This indicates spectra from that source were generally distinct from other spectra included in the initial spectral library.

When spectral libraries were used to model cover fractions across the Rim Fire, spectra generated from the AVIRIS images were disproportionally selected, given their abundance in the reduced spectral libraries (

Table 5. Percentage of pixels modeled by source for each endmember reduction technique. A single pixel can be modeled by multiple different sources, resulting in totals that exceed 100%.

Spectra Source	EMC	uSZU EMC	In-CoB	uSZU In-CoB	IES	uSZU IES	Reduced IES	uSZU Reduced IES
Not Modeled	13%	18.3%	16.8%	20.9%	7%	8%	14.3%	13.6%
AVIRIS	87%	81.7%	42.1%	53.2%	61.8%	55.3%	48.6%	32.5%
ASD	0	0	15.7%	6.2%	5.4%	10.3%	6.9%	26.4%
Wind River	0	0	45.8%	38%	47.1%	44.4%	47.6%	41.6%

). All spectral sources were used to some degree. It should be noted that not all endmember classes modeled in this study were collected for each spectral source and there were large differences in the number of spectra collected for each source.

3.5. Validation

3.5.1. WorldView-2 Based Validation

The classified WorldView-2 imagery and the AVIRIS-based MESMA cover fractions were compared using the 120 randomly selected polygons across the two image types (

Table 6. The coefficient of determination (r²), intercept, and slope value of a linear fit between WorldView-2 estimated cover fractions as the dependent variable and AVIRIS derived MESMA cover fraction for different endmember and band selection techniques as the independent variable. Green vegetation is abbreviated as GV, non-photosynthetic vegetation as NPV.

		EMC	uSZU EMC	In-CoB	uSZU In-CoB	IES	uSZU IES	Reduced IES	uSZU Reduced IES
Char	r2	0.605	0.741	0.727	0.642	0.620	0.538	0.687	0.594
	Intercept	0.355	0.303	0.212	0.202	0.321	0.223	0.255	0.310
	Slope	0.721	0.755	0.841	0.835	0.765	0.783	0.771	0.750
GV	r2	0.750	0.770	0.836	0.871	0.848	0.846	0.853	0.861
	Intercept	−0.054	−0.023	−0.026	−0.065	0.020	0.000	-0.010	−0.019
	Slope	0.653	0.769	0.708	0.675	0.895	0.812	0.804	0.748
NPV	r2	0.086	0.099	0.164	0.249	0.209	0.237	0.165	0.245
	Intercept	0.110	0.103	0.109	0.126	0.099	0.107	0.099	0.094
	Slope	0.240	0.521	0.426	0.375	0.968	0.943	0.377	0.699
Soil	r2	0.273	0.261	0.088	0.042	0.014	0.049	0.075	0.057
	Intercept	0.013	0.015	0.015	0.013	0.019	0.016	0.011	0.012
	Slope	0.641	0.102	0.222	0.965	0.041	0.088	0.2	0.103

). The best linear correlations between WorldView-2 fractions and MESMA were observed for the GV cover, which displayed a near 1:1 relationship and had generally high r² values regardless of the endmember selection technique used. Correlation coefficients for char were generally fairly high (between 0.59 and 0.741), however, all models showed a tendency to under-model char, as demonstrated by linear models consistently having an intercept above 0.2. Correlations were generally poor for NPV and soil, with r² values often below 0.2. Larger endmember libraries did not always translate to higher fractional cover accuracy; for example, the highest r² for the char fraction was found for the uSZU EMC library (0.741).

Plots of the linear regression for the uSZU IES MESMA run and the WorldView-2 classifications are provided below (Figure 6). The plots are broadly representative of the trend of most MESMA versus WorldView-2 relationships. There is general support for the linear relationship of the GV identification by MESMA and the WorldView-2, with a high r² value and little systematic error. Based on the manual interpretation of the WorldView-2 imagery, MESMA appears to be systematically modeling a lower fraction of NPV and char cover, and modeling a higher fraction of soil cover.

3.5.2. Comparison with GeoCBI Plot Data

Due to the similarity between MESMA unmixed models used in this study, we limit our comparison to only the cover fractions generated from MESMA using the IES library and the full AVIRIS bands. First, the relationship between individual cover fractions and GeoCBI was assessed (Figure 7). The relationship between the GV fraction and GeoCBI appeared to be inversely linear (a linear regression produced an intercept of 2.93, a slope of −1.89 and r² value of 0.644). In contrast, the relationship between NPV and GeoCBI is clearly non-linear. NPV fractions are generally highest at GeoCBI values of 1.5 to 2.5 suggesting that NPV is high at moderate severities, but low at the highest and lower severity levels. The char fraction appears to be more of a binary relationship only appearing in the model at the highest GeoCBI levels (above 2.75) and is at or near zero for all other GeoCBI levels.

3.5.3. Comparison with Field Tree Status Data

We also compared the MESMA cover fractions of GV, NPV, and char with the percentage of canopy level trees with mostly green, brown, and black needles (Figure 8). Similar to the GeoCBI values, the char percentage and percentage of black trees appeared to have a near binary relationship in the plots. If the percentage of field identified black trees within the plot was under 50%, then typically no char fraction was modeled for the pixel, if it was above 50%, then the char fraction was almost always modeled as close to 100%. The relationship between the brown trees and NPV appeared to generally be more linear, although the modeled NPV cover was never greater than 50% even with near 100% brown trees. The relationship between the green trees and green vegetation may be positively linear, but the amount of scattering makes interpretation difficult. It is important to remember that few AVIRIS pixels are purely tree cover, most are a mixture of tree and substrate, so a pixel with 100% brown trees may indeed be made up of 50% other materials. Finally, it is important to note that the percent of trees in the plot may not be an ideal reference for the cover fraction as it does not account for variation in tree density, and this may partially contribute to the observed biases.

4. Discussion

4.1. Potential Bias and Uncertainty in the Cover Types

Whereas previous studies found the relationship between the char fraction and GeoCBI to be linear [18,58] this was not the case in this study. The GV fraction retrieval by MESMA, was found to be linear, however. GV’s high accuracy is likely due to both the spectral separability and the temporal persistence of the class. In contrast to GV, at the spectral resolution of WorldView-2 soil, NPV, and ash were all somewhat difficult to separate. Human interpretation likely reduced inter-class confusion in classifying the WorldView-2 imagery, but this was still a source of error associated with the WorldView-2 image classification. In particular, separation of the soil and ash areas was particularly difficult with the WorldView-2 imagery. This is partially because of insufficient spectral information and partially because even at 2 m pixel resolutions, many pixels were likely a mixture of mostly ash with some soil. This mixture would impact the spectral signal, and therefore be reflected in the MESMA generated cover fractions, but these areas were generally just classified as ash in the manual interpretation of the WorldView-2 image. Another possible factor in diminished char is the amount of time between when the fire burned the area of the GeoCBI plots, which was in late August and early September (evaluated from fire spread map in [59]), and image acquisition, which was 28 September for the WorldView images and 17 November for AVIRIS. During the time period between the fire occurring and AVIRIS image acquisition, there were four rain events totaling 64.8 mm at a weather station in Groveland, CA (data accessed from http://www.ncdc.noaa.gov/cdo-web/datasets#GHCNDMS) within 20 km from where the field plots were collected. Combined with wind, rain likely attenuated the ash portion of the char signal. The impact of weather, combined with other previously mentioned factors, likely account for why char was more prevalent and soil less prevalent in the WorldView-2 classification than the MESMA cover fraction (Figure 6).

4.2. Evaluation of MESMA Techniques’ Performance

The balance of library complexity with accuracy must be considered if MESMA is to be used on spaceborne imaging spectroscopy at regional to global scales for ecological monitoring [58]. Schaff et al. [31] and Roth et al. [32] have shown IES endmember selection produces high classification accuracy in the two endmember MESMA case. However, IES selected the most endmembers in its final spectral library in our study and when compared to other methods in [32]. In this study, regardless of the technique used to generate cover fractions, the more spectrally unique classes (GV and ash) had approximately similar performances when compared with the WorldView-2 cover fractions. Using uSZU In-CoB to generate MESMA cover fractions produced relatively high r² values when compared to WorldView-2, even though the number of endmembers used and processing times were significantly less. That a variant of In-CoB would perform as well as IES is not consistent with [32]. One possible explanation for this is that when using endmembers for classification of vegetation types, as was done in [32], some degrees of mixture within endmembers, particularly with soil and NPV, is necessary for accurate classification across a landscape. If endmembers within the endmember classes of the starting library are to some degree mixed, in order maximize kappa IES will tend to select mixed pixels. However, when using MESMA for spectral unmixing to estimate cover fractions, as was done in this study, the purest endmembers will obtain the highest accuracy.

The margin of computational efficiency created by uSZU, a 30–50% reduction in processing times, was similar to Somers and Asner [47]. However, unlike [47], we did not observe a clear trend of increased accuracy (Table 6). In Figure 5 and Table 3, it is clear that the techniques that used uSZU generally mapped less NPV than techniques that did not. Spectral regions that are critical for the separation of soils and NPV due to their association with lignin and cellulose absorption, such as the 2300–2400 nm region [60], did not get selected by uSZU (Figure 4). Many of the bands which were selected appear associated with atmospheric noise. In Veraverbeke et al. [58] bands were reduced to similar numbers as the uSZU, however they were not selected in a way that was designed to maximize differences between endmember classes, resulting in a much larger decrease in accuracy than was observed here. Future studies should investigate other band selection techniques and the optimal bands for the simultaneous separation of GV, NPV, soils, and ash.

4.3. Endmember Sources and Endmember Selection

Beyond pre-processing, another challenge to SMA becoming a global means of objective comparison is determining a proper base spectral library. Studies have shown the timing of the acquisition of the base spectra and the image is important [61], as is the spatial scale that the spectra were acquired [35,51]. It is notable that spectra collected in a different spatial location than the images and at different spatial scales than the images, were frequently selected to model cover fractions. This suggests that a common spectral library could be developed to map fire across at least a regional scale. It also suggests that even in environments not composed of as many complex materials as the urban [35,62,63], a diversity of spatial scales may be beneficial to cover fraction mapping.

Although, impurity of endmembers collected in non-laboratory settings is inevitable, it is a particular challenge at the spatial scale of AVIRIS or proposed spaceborne imaging spectrometers. It is notable that the endmember classes which it is most difficult to get pure pixels from images, NPV and soil (that was not exposed rock), also appeared least accurate. Although spectra collected in the field were included in our spectral library, these may not perfectly scale to the canopy level AVIRIS observations [51]. AVIRIS-Next Generation (AVIRIS-NG; [64]) poses a potential means for overcoming some of these challenges [55]. The superior spatial and spectral resolutions allows for collecting and evaluating endmembers that are closer to being pure while still being collected at the canopy-scale of AVIRIS.

4.4. SMA as a Novel Means for Assessing Fire Severity

Since properly identified cover fractions have physical meaning, they have the potential to be an objective and global means of assessing fire severity. Cover fractions, if shown to be accurate and comparable to field derived fire severity estimates, have the potential to overcome many of the criticisms of NBR-based indices [10,19,65,66]. With the potential launch of several spaceborne imaging spectrometers, and considering the demonstrated higher correlations of SMA with imaging spectroscopy data compared to broad-band data [58], the use of SMA for fire may increase in the future. Since the variation of soil and GV’s spectral profile across space is one of the key reasons for NBR’s subjectivity, if SMA is going to become an important means of assessing burn severity, endmember variability must be accounted for, and MESMA is one of the most reliable technique for doing so [22]. In this study, only a single fire is included in the study area; however, if SMA is to become a commonly used tool for fire severity assessment, a global assortment of geographic locations and their spectral variability will need to be tested and compared. As imaging spectroscopy data becomes more common, both through more airborne and spaceborne acquisitions, global assessment will become viable.

One possibility for assessing change using SMA is using a differenced SMA (dSMA) approach. By using dSMA the entirety of the information provided by MESMA classes could be used, potentially providing a robust and more ecologically meaningful method of evaluating fire severity. However, this approach would have an inherent disadvantage of all differenced imagery comparisons, in that relatively analogous pre-fire imagery would be needed. Since all imaging spectroscopy data with signal-to-noise ratios and spectral and spatial resolutions similar to AVIRIS are currently acquired from airborne platforms, pre-fire data is rare; however proposed spaceborne imaging spectrometers such as HyspIRI and ever-increasing computational power would make this analysis possible at regional to global scales. Future studies evaluating dSMA usefulness in fire severity evaluation will require careful planning in field validation plot placement to assure that plots are placed in areas that pre-fire were heterogeneous in terms of type and percentages of vegetation cover. Given the relative accuracy of mapping the GV cover class, using a differenced GV cover fraction (dGV) is likely to have high correlation with the field measures of burn severity. In addition, the high correlation observed for the char cover type, suggests that adding a post-fire char cover to a dGV assessment could further enhance the discrimination of burn severity.

5. Conclusions

We demonstrated the utility of imaging spectroscopy combined with MESMA for fire severity mapping over a large fire in California’s Sierra Nevada. While currently limited to airborne acquisitions, future spaceborne missions will allow large-scale application of these techniques. One aspect that will be important for processing imaging spectroscopy and MESMA globally are techniques that maximize accuracy while minimizing calculation time. In this study, we performed MESMA using a spectral library that included several different collection methods and multiple different spatial scales. We demonstrated that band reduction can significantly reduce computational time with only small differences in performance. We also evaluated the performance of several endmember selection techniques and found that these can also be optimized between performance and calculation time. Further research is needed that evaluates cover fractions relation to fire severity comparing multiple fires across regions. However, the identification cover fractions represents a potentially objective and physically meaningful evaluation of fire severity using remote sensing.