Endmember bundles - IEEE Xplore

IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 38, NO. 2, MARCH 2000

1083

Endmember Bundles: A New Approach to Incorporating Endmember Variability into Spectral Mixture Analysis C. Ann Bateson, Gregory P. Asner, and Carol A. Wessman

Abstract—Accuracy of vegeation cover fractions, computed with spectral mixture analysis, may be compromised by variation in canopy structure and biochemistry when a single endmember spectrum represents top-of-canopy reflectance. In this article, endmember variability is incorporated into mixture analysis by representing each endmember by a set or bundle of spectra, each of which could reasonably be the reflectance of an instance of the endmember. Endmember bundles are constructed from the data itself by an extension to a previously described method of manually deriving endmembers from remotely sensed data. Applied to remotely sensed images, bundle unmixing produces maximum and minimum fraction images bounding the correct cover fractions and specifying error due to endmember variability. In this article, endmember bundles and bounding fraction images were created for an airborne visible/infrared imaging spectrometer (AVIRIS) subscene simulated with a canopy radiative transfer/geometric-optical model. Variation in endmember reflectance was achieved using ranges of parameter values including leaf area index (LAI) and tissue optical properties observed in a North Texas savanna. The subscene’s spatial pattern was based on a 1992 Landsat Thematic Mapper image of the study region. Bounding fraction images bracketed the cover fractions of the simulated data for 98% of the pixels for soil, 97% for senescent grass, and 93% for trees. Averages of bounding images estimated fractional coverage used in the simulation with an average error of 0 05, a significant improvement over previous methods with important implications for regional-scale research on vegetation extent and dynamics. Index Terms—Canopy radiative transfer, Dikin’s affine algorithm, endmember, endmember bundle, endmember variability, simulated annealing algorithm, spectral mixture analysis.

I. INTRODUCTION

F

RACTIONS of ground cover computed from a spectral unmixing of multispectral and hyperspectral datasets have been critical for accessing ecologically significant plant attributes and processes hidden in remotely sensed data. Adams et al. [1] classified four Landsat thematic mapper (TM) images of the Brazilian Amazon using spectral mixture analysis, and Manuscript received September 24, 1998; revised April 20, 1999. This work was supported by NASA Innovative Research Grant NAGW-4689 and NASA Land Use/Land Cover Change Grant NAG5-6134. C. A. Bateson is with the Center for the Study of Earth from Space (CIRES), University of Colorado, Boulder, CO 80309-0216 USA (e-mail: [email protected]). G. P. Asner is with the Department of Geological Sciences and Environmental Studies Program, University of Colorado, Boulder, CO 80309-0216 USA. C. A. Wessman is with the Department of Environmental, Population and Organizm Biology and is also with the Center for the Study of Earth from Space (CIRES), University of Colorado, Boulder, CO 80309-0216 USA. Publisher Item Identifier S 0196-2892(00)02484-0.

revealed changes in land use related to such processes as forest clearing and regrowth. A spectral unmixing of airborne visible/infrared imaging spectrometer (AVIRIS) imagery produced cover fractions discerning management practices of burning and grazing in a Kansas tallgrass prairie [2]. Recently, spectral mixture analysis of a TM scene acquired over North Texas was combined with a one-dimensional (1-D) canopy radiative transfer (RT) model inversion to simulate a three-dimensional (3-D) RT model retrieval of leaf area index (LAI) from multiangle advanced very high resolution radiometer (AVHRR) data [3]. These LAI estimates were instrumental in computing ecologically important vegetation parameters of total diurnal fraction of absorbed photosynthetically active radiation (FAPAR) and daily carbon uptake. Extensive and intensified land use is causing significant changes in the structure and function of terrestrial ecosystems with consequent effects on biosphere–atmosphere interchanges. Remote sensing is required for large scale measurements needed to monitor structural change over time and to analyze and assess functional effects via ecosystem models [4]. The above applications of spectral mixture analysis have revealed its potential for providing a suite of spatially distributed variables that can initialize and constrain these ecosystem models. Since exact rather than relative measures of cover fractions (e.g., NDVI) are mandated by the models and by structural analyses of land use/cover change, methods to either significantly increase or gauge the accuracy of the fraction estimates become a necessity. Traditionally, spectral mixture analysis models a reflectance spectrum as the linear combination of a finite number of spectrally unique signatures of pure ground components, referred to as endmembers. With an ideal choice of endmembers, the coefficients in the linear combination are nonnegative, sum to 1, and therefore, interpretable as cover fractions. Such a set of ideal endmembers together with the cover fractions is called a mixture model for the spectrum. A mixture model based on three or four endmembers has the simple geometrical interpretation as the triangle or tetrahedron whose vertices are the endmembers. Cover fractions are determined by the position of the modeled spectrum within the interior or on the boundary of the figure. Noninstrumentally induced errors can be injected into cover fractions because multiple scattering of photons between different ground components invalidates the linearity assumption of the model [5], [6]. One strategy to ameliorate the impact of nonlinear mixing on spectral mixture analysis attempts to linearize reflectance by conversion to single scattering albedo [7],

0196–2892/00$10.00 © 2000 IEEE

1084


Fig. 1. Scatterplot of 2-D spectral data. The plot illustrates the effects of endmember variability on fitting a triangle to the data. Points corresponding to spectra representative of trees, grass, or soil are contained within the labeled regions. The vertices T, G, and S of the smallest triangle containing the data are not within these regions and hence do not correspond to the spectral reflectance of trees, grass, and soil. In fact, they may represent unrealistic spectra with negative or superunity reflectance and unmixing with the three vertices T, G, and S, which will introduce error into the cover fractions.

[8]. In [9], a Monte Carlo ray trace model found that nonlinear mixing, which occurs when a light ray interacts with multiple elements, may be minimized by acquiring images from view angles close to the hot spot direction (that is, from angles where a majority of the scattered photons are not seen). However, probably the most profound source of error in linear spectral mixture analysis lies in biophysical and biochemical variability of vegetation types leading to concomitant variability in top-of-canopy reflectance [10], [11]. Fig. 1 illustrates how fractional errors can occur under the assumption that all pixels are mixtures of three endmember spectra (soil, grass, and tree) when in fact, there are spectral variations in the endmembers. Under the above assumption, the endmembers should be the vertices of a triangle encompassing their mixtures, but there is no triangle encompassing all the mixtures in Fig. 1 with vertices in regions representing the endmembers. In order to incorporate endmember variability into spectral mixture analysis, we represent a landscape component type not with one endmember spectrum but with a set or bundle of spectra, each of which is feasible as the spectrum of an instance of the component (e.g., in the case of a tree component, each spectrum could reasonably be the spectral reflectance of a tree canopy). A pixel may have infinitely many mixture models with endmembers selected from the bundles. Bundle unmixing finds the maximum and minimum endmember fractions in those models and hence the range of endmember fraction values possible for a pixel. Previous efforts to acknowledge endmember variability have been based on image-wide estimates of uncertainty [10], fuzzy membership of a mixed pixel in various possible mixture proportions [12], foreground-background analysis [13], and unmixing with multiple field or image endmembers [14], [15],[16]. Unlike methods to compute image-wide error estimates, the fuzzy fraction approach and the endmember bundle method recognize that uncertainties in the fractions

vary from pixel to pixel. For example, in Fig. 1, there is no uncertainty in the fractions of a mixture on the boundary of the triangle, while for mixtures in the center, there are many mixture models based on endmembers selected from the tree, grass, and soil bundles. The fuzzy fraction approach uses image endmembers and produces information (fuzzy memberships) supplemental to the fraction range determined from a bundle unmixing. The foreground-background method determines the projection in the hyperspace of the spectral data that minimizes variation both in the endmember of interest (foreground) and all other endmembers (background). Therefore, it works best when almost parallel hyperplanes can be well fitted to the foreground and background spectra, respectively. Otherwise, significant error will persist. Although the method of multiple endmember spectral analysis could be applied to a collection of image endmembers, currently, the technique uses a library of endmembers collected from the field and pixel by pixel tries all combinations of endmembers in search of the best fit. This method is very intensive in terms of both human and processor time. Finally, the method of Bajjouk et al. [17] is similar to the method of this article if a proliferation of endmembers in a 2-D or 3-D space is viewed as really only several endmembers with variability (i.e., bundles). However, the technique we use to derive endmembers and bundles does not limit us to a 2-D or 3-D space and is more systematic than Bajjouk et al. [17]. In [18], we used a nonlinear optimization algorithm to determine the range of possible fraction values of a ground cover type in a pixel. In the current article, we use essentially the linear programming method of Bajjouk et al. [17] to more efficiently compute maximum and minimum fraction images. Also, while endmember spectra have been collected from pure pixels in the image or from the field to estimate the spectral variability of a vegetation type, a salient feature of our method is the construction of the bundles from the image itself as an extension to a previously described method of manually deriving endmembers from spectral data [19]. Image-derived endmembers have an advantage over image endmembers in that they represent pure endmember spectra when there are perhaps no homogeneous pixels in the image. Moreover, unlike field endmembers, derived endmembers obviate the need to do extensive ground work to develop a spectral library. Endmember bundles, which are developed from a set of derived endmembers, inherit these advantages. The main goal of this article is to present a method for constructing endmember bundles from spectral data and performing bundle unmixing. This method was tested by constructing endmember bundles and bounding fraction images for an AVIRIS subscene simulated with a hybrid plant canopy radiative transfer/geometric-optical model. In several studies, this model has been used to produce realistic spectra from realistic values of field-measured vegetation properties (Table I) [3], [11], [20]. Simulated rather than remotely sensed data were used for testing the methods in order to have completely accurate “ground truthing” for a large number (40 000) of pixels and to avoid the need for image calibration or correction. With this data set, the bounding fraction images were tested to see how well they bracketed the cover fractions used in the simulation. For each cover type, minimum and maximum bounding images

BATESON et al.: ENDMEMBER BUNDLES

TABLE I RANGE OF PARAMETER VALUES FROM A NORTH TEXAS SAVANNA USED TO SIMULATE VARIABILITY WITH CANOPY RADIATIVE TRANSFER/GEOMETRICAL-OPTICAL MODEL

were averaged to produce estimates of cover fractions. We also present a method of perturbing these averages in order to obtain fractions lying between the bounding images and summing to 1. In addition, for comparison, we performed a traditional unmixing with image endmembers, derived endmembers, and the mean of endmember spectra created by the canopy RT model for the simulation. For all methods of selecting endmembers and for the bundle means and perturbed means, comparison was made between computed fractions and the cover fractions used in the simulation. However, when applied to unsimulated data, only bundle unmixing provides per-pixel estimates of the errors in its computed fractions. II. METHODS A. Model Simulation A canopy radiative transfer/geometric-optical model was used to simulate an AVIRIS image based on the variability of plant tissues, canopy architecture, and landscape structure measured in a heterogeneous North Texas savanna ecosystem (Table I), [3]. Details of the model have been provided by Iaquinta and Pinty [21], Li and Strahler [22], Asner et al. [23], [3], and Asner and Wessman [24], but a brief summary follows. Individual pixels are populated by tree and senescent grass canopies with some fraction of bare soil. Scattering of incident

1085

photons by each tree and grass canopy is modeled using a multistream approach by which photons are permitted to travel in 96 directions, or 12 directions per octant of a unit sphere. Both leaf and nonphotosynthetic vegetation are included in the radiative transfer formulation [24], [25]. Single scattering of photons by the mixed canopies is solved exactly, while multiple scattering is simplified to a single-angle (sun zenith) problem [21]. The model is designed specifically for simulating hyperspectral data, as wavelength-independent calculations (e.g., LAI) are made only once per simulation, while those calculations requiring tissue optical properties (e.g., multiple scattering) are iterated by wavelength. The model produces pixel-level reflectance values from the following parameters: tree and grass LAI, tree and grass leaf angle distributions (LAD), leaf optical properties, soil reflectance, stand density, tree crown dimensions, sun and view geometry, and a hotspot parameter for each vegetation type. LAI is given on an m m basis, LAD can be modeled as erectophile, planophile, plagiophile, and uniform [26]. In this study, LAI was randomly selected within a range measured in the North Texas savanna (Table I), [3]. Leaf angle distributions were modeled as ellipsoidal with random selection of the mean inclination angle. Scattering characteristics at the tissue and soil level were modeled as isotropic, and here, were parameterized using measurements obtained in the savanna region in 1996 [3]. Both hemispherical leaf reflectance and transmittance properties of plants were used in the radiative transfer simulation. The tissue spectra were collected from the Texas savanna region as part of a full spectral-range (400–2500 nm) database [11]. The optical properties of both sides of leaves were accounted for in the simulation. For all analyses in this study, solar zenith and azimuth angles were set, respectively. View zenith and azimuth angles were both set to 0. The four endmember groups were tree, senescent grass, shade, and bare soil. Green grass was not used in the simulation because grass in the study site at the time of the TM acquisition was mostly senescent. Variability within each endmember was a function of the tissue optical and canopy structural characteristics already mentioned (Table I). The normalized difference vegetation index (NDVI) from a 1992 Landsat TM image of the Texas savanna region was used to create a nonrandom spatial organization in the simulated AVIRIS scene. Tree stand density values were scaled to the highest NDVI value in the (ndvi / )NDVI) + Landsat scene [stand density 0.30]. A crown geometric-optical (G-O) submodel computed tree, shade, and sunlit understory cover fractions using this stand density value, randomly selected crown dimensions, and prescribed sun-view geometry [22]. Soil and senescent grass fractions were then derived from the G-O sunlit understory NDVI , the senescent grass fraction fraction. For was scaled to the minimum NDVI [grass = sunlit fraction ndvi NDVI-0.4)], and the sunlit soil fraction was the residual value needed to sum the four end, the sunlit grass fraction was members to 1. When NDVI calculated as a random fraction of the G-O sunlit understory fraction [grass = sunlit fraction * random value (0.0–1.0)]. Noise was added to each spectrum based on the published AVIRIS in-flight SNR over Ivanpah Playa, CA, on the June 4,

1086


1996 [27]. To compute the noise component of a spectrum, the reflectance in each band was divided by the SNR for that band and a mean dark-signal offset was added. The fraction images from the AVIRIS simulation model will be referred to in the sequel as the ground truth (GT) fraction images with one exception. The high correlation between stand density and shade fractions prevented separation of trees and shade. Thus, their sum was considered as a single endmember in the analysis (called the “GT tree fraction”). B. Linear Spectral Mixture Analysis Linear spectral mixture analysis models the reflectance spectrum of each pixel in a multispectral/hyperspectral image in terms of endmember reflectance according to the equations and constraints below (1)

(2) is the fractional coverage of endmember in the where pixel, ranges over the spectral bands of the instrument, and is the residual error in band . If the number of endmem, then in (1) and (2) above, bers equals the number of bands there are equations in unknown, which can be uniquely inverted to solve for the ’s with no error terms in the equations and fractions are a of (1). In this case, the endmembers mixture model for the pixel, provided that the ’s are nonnegative. To derive endmembers and endmember bundles from the data, we first perform a principal component (PC) analysis on the mean-corrected1 data to reduce its dimension to the number of components accounting for most (%) of the spectral variance [19]. In this PC-reduced space, we select endmembers, construct bundles, and perform bundle unmixing. C. Bundle Construction The schematic diagram in Fig. 2 summarizes the sequence of processing steps described in this section to construct endmember bundles. After the dimension of the data has been reduced through a PC analysis, the next step in bundle construction is to use an automated method to find a simplex2 in the PC-reduced space enclosing all data points and satisfying a specified criterion. One possible criterion requires searching for the minimum volume simplex enclosing the data [28]. For this paper, we used the simulated annealing algorithm in the Appendix to find the minimum volume simplex constructible from the data cloud by a procedure described in the Appendix. If the spectral data had two dimensions as in Fig. 1, then the endmembers found with the simulated annealing algorithm would be the vertices of a triangle. -dimensional space As explained in [19], points in the determined by the PC analysis can also be viewed as spectra in 1That is, all spectra are averaged, and this average is subtracted from each spectrum in the image 2A simplex is a triangle in the plane, a tetrahedron in 3-D space, and in -space, the analogous geometrical figure that has vertices and consists of all mixtures of those vertices.

N

N +1

Fig. 2. Schematic representation of bundle construction and unmixing. The AVIRIS data were projected into the space determined by the principal components accounting for most (> 99%) of the variance. A simulated annealing algorithm derived endmembers (Fig. 3) as vertices of a simplex, enclosing the data in the PC-reduced space. The manual endmember selection method relocated these endmembers (Fig. 4) so that their associated spectra were acceptable as reflectance spectra of ground components. The relocated spectra were used as seeds to grow the endmember bundles (Fig. 5) by including in each bundle realistic spectra sufficiently correlated to its seed. Maximum and Minimum fraction images are the products of a bundle unmixing.

the full set of spectral bands. Because of endmember variability (Fig. 1) and nonlinear mixing, automatically generated spectra may not be acceptable as the reflectances of ground components. In fact, they may have negative or superunity reflectances. The manual endmember selection method [19] initialized with automatically generated endmembers, can be used to find spectra acceptable as endmembers. With reference to Fig. 1, this step moves our candidate tree, grass, and soil endmembers , and to points within the from their respective vertices regions labeled tree, grass, and soil. Once acceptable endmembers are found, they are used in the PC-reduced space as seeds to grow bundles in the following way. Initially, each bundle contains its seed. Nearby points are included in the bundle as long as their corresponding reflectance spectra have values ) between 0 and 1 and are sufficiently correlated ( to the spectrum associated with the seed endmember. As an option, correlation coefficients may be computed over the entire spectrum or over subranges of reflectance values individually. Other options include specifying reflectance limits for each bundle (e.g., tree bundle spectra should never have reflectance in any spectral band greater than 0.6) and the Euclidean distance between points in the PC-reduced space consecutively tested for similarity of their associated spectra to the seed endmember’s spectrum. Finally, the bundle is the


1087

convex closure of all points tested in the PC-reduced space whose associated spectra satisfy all reflectance limits and correlation requirements. Equating the bundle with the convex closure means that the linear combination of any two points in the bundle with nonnegative coefficients summing to 1 (i.e., any mixture of the two points viewed as endmembers) also belongs to the bundle. That is, we are recognizing that there may be endmember variability within a pixel as well as from pixel to pixel. In the bundle construction, there is a minimum finite set of points, called the convex hull, whose mixtures determine the bundle.

TABLE II ENDMEMBER SELECTION AND UNMIXING METHODS FOR THREE TRADITIONAL MIXTURE ANALYSES OF THE SIMULATED DATA. ALSO SPECIFIED ARE THE TREATMENTS OF SHADE FRACTIONS IN ESTIMATES OF GROUND TRUTH TREE SHADE COVERAGE

+

D. Bundle Unmixing Below are the bundle unmixing equations and constraints that model a spectrum’s projection into the space of the first principal components where each . We replace (4) with (3)

(8)

(4)

is a very large number forcing to equal 0 when where is initialized to and (8) is minimized [30]. If each is initialized to 1, then the equations in (7) are true. When (8) is minimized, the equations in (7) will coincide with those in (3) since s will equal 0. Also, since linear programming always finds nonnegative values for the variables, the inequality in (4) is computed by miniwill hold automatically. Hence, is similarly mizing (8) subject to the constraints in (7). substituted for computed by minimizing (8) above with and subtracting the minimum from 1. In addition to computing maximum and minimum fraction images through bundle unmixing, the mean fraction image for each bundle is also computed as

where each bh belongs to the convex hull of a bundle, and is the total number of points in all bundle hulls. In the presis greater than , and (3) ence of endmember variability, and (4) have infinitely many solutions. Bundle unmixing uses a slight variation of the linear programming method of Bajjouk et al. [17], which is based on Dikin’s affine algorithm method and minimum [29], to compute the maximum of all possible fractions of bundle in a pixel. For example, if is the hull of the tree bundle, then, to find for a pixel, bundle unmixing would minimize the sum below of all tree hull fractions subject to the linear constraints in (3) and (4) (5) where (6) in (5) has the effect of removing from the sum The term those fractions not relating to the tree bundle because they are , each is replaced by multiplied by 0. To find in (5) above, so that now only fractions not relating to the tree bundle are left and their sum is minimized subject to the constraints in (3) and (4). This minimum is subtracted from 1 to . That is, the tree fraction is maximized in a pixel obtain by minimizing the area not covered by trees. Dikin’s affine linear programming algorithm is an iterative algorithm that needs to be initialized with a set of ’s satisfying the equations and constraints in (3) and (4). To facilitate initialization of the algorithm, we introduce the dummy variable and replace the first set of equations in (3) with

(7)

(9) For applications requiring fractions summing to 1, while still lying within the limits of the bounding fraction images, we suggest the following computation. For each pixel, the sum of is computed. The mean its mean fractions and sum of half intervals are used to for bundle compute a perturbed fraction (10) The images formed from the perturbed fractions are called perturbed mean bundle images. E. Traditional Unmixing For comparison, the simulated AVIRIS image was unmixed with image endmembers and photometric shade3 using the weighted least squares (WLS) spectral mixture analysis algorithm in [31] to minimize the sum of squares of the errors in (1) subject to the constraints in (2). In this unmixing, each image endmember was the spectrum of the image pixel with the highest fraction in the GT image for the given cover type. Photometric shade was used as a surrogate for tree shade, since the pixel with the highest shade fraction was also the pixel 3Reflectance

in all bands is 0.

1088


Fig. 3. Five endmembers derived from the simulated AVIRIS data set using an algorithm based on simulated annealing. This algorithm estimated the vertices of the minimum volume simplex containing the data in a four-dimensional (4-D) PC-reduced space. That is, the endmembers above are the 4-D analog of the vertices T ; G, and S in Fig. 1.

with the highest tree fraction. The simulated image was also unmixed by the method in [19], with endmembers derived with the simulated annealing algorithm. Because these derived endmembers contain the spectral data in the simplex, they determine in the PC-reduced space, all pixel fractions are nonnegative and sum to 1. Finally, the simulated image was unmixed with mean endmember spectra created with the plant canopy RT model. More specifically, endmember spectra created with the plant canopy RT model for the simulation were saved, and spectra for each ground cover type were averaged to obtain a set of four endmembers (tree, senescent grass, shade, and soil). When the simulated data were unmixed with this endmember set, the senescent grass fractions were all negative. However, better results were produced by reducing the endmembers to the space of the first four principal components of the simulated data. Again, unmixing was performed with the WLS method in [31] to insure that fractions satisfied the constraints in (2). The means (in PC-reduced space) used in this unmixing will be referred to as GT mean endmembers, since they are analogous to the mean of a set of field spectra totally representing the variability of a given vegetation or ground cover type. Table II summarizes the methods of endmember selection and spectral mixture analysis considered in this article.

III. RESULTS AND DISCUSSION A PC analysis of the data found that the first four components explained over 99% of the variance and that the simulated data were mean-corrected and projected into the space determined

by those components. In this PC-reduced space, the five endmembers in Fig. 3 (one tree, two senescent grass, a shade, and soil endmembers) were derived with the simulated annealing algorithm. The manual endmember selection tool relocated the endmembers in Fig. 3 to obtain the endmembers in Fig. 4. Note that the tree endmember no longer has negative reflectance, the shade endmember resembles a vegetation endmember, and the soil endmember has a smoother appearance in the 0.5–1.0 m range. These endmembers were used as seeds to grow the endmember bundles. The two senescent grass bundles were constructed by looking at correlations computed on the three uninterrupted wavelength intervals separated from each other by the water bands,4 specifically, the intervals 0.4–1.26, 1.43–1.76, and 2.0–2.5 m. Since water bands do not give information on surface reflectance, they have been removed from all spectra, a customary practice with AVIRIS data [32]. For the tree and shade bundle, we separated the visible and NIR wavelengths and divided the spectral bands into the four wavelength intervals 0.4–0.7, 0.71–1.26, 1.43–1.76, and 2.0–25 m. Correlations with the seed tree endmember were computed separately on the four intervals. For both grass and tree bundles, the intervals were selected to separate spectral features from each other in order to help prevent a very good fit in one part of the spectrum from compensating for a poor fit in another spectral region. The maximum reflectance values of the tree, senescent grass1, and senescent grass2 spectra were restricted to be no more than 0.60, 0.80, and 0.65, respectively. These reflectance limits were based on observed field values [11], [33] and the appearance of unnatural features in spectra whose reflectance exceeded the 4That is, the AVIRIS bands were dominated by atmospheric water absorbtion.


1089

Fig. 4. Five endmembers found by relocating the endmembers in Fig. 3 with the manual endmember selection tool in order to find spectra within the endmember bundles.

Fig. 5.

Spectra on the convex hull of the five bundles grown from the endmembers in Fig. 4.

maximum. Also, we restricted the soil bundle to consist of a single spectrum under the assumption that there was only one soil type. Some spectra were removed from the bundles on the

basis of a visual inspection because they appeared as mixtures or they had unnatural features. The final endmember bundles are shown in Fig. 5.

1090


Fig. 6. Minimum and maximum bounding images bracketing the ground truth images, which contain the fractions determined in the simulation of the AVIRIS data set.

Because of the similarity between tree and shade endmembers, we combined their fractions in the bundle unmixing as well as the fractions of the two senescent grass endmembers. Two and may be combined in a bundle endmember bundles unmixing by redefining in (6) by

TABLE III THE DIFFERENCE BETWEEN MAXIMUM AND MINIMUM BOUNDING IMAGES AVERAGED OVER ALL PIXELS AND FOR EACH ENDMEMBER. STANDARD DEVIATIONS ARE ALSO SHOWN

(11) Resulting maximum and minimum images are shown in Fig. 6, bracketing the GT images. These maximum and minimum images bracket the ground truth fractions for 98% of the pixels for soil, 97% for senescent grass, and 93% for tree fractions values. That is, the bundle construction method captured almost all of the endmember variability. The minimum and maximum images define for each pixel an interval whose length provides an estimate of the error in the fractions when wrong endmembers are used to unmix the pixel. Table III gives the length of these intervals averaged over all pixels and for each ground cover type. If we estimate each endmember fraction in a pixel with its value in the mean bundle image, then on average, the estimated error in tree fractions , senescent grass is is , and soil is . Figs. 7 and 8 show, respectively, the three endmembers selected from the image and the four mean GT endmembers. Photometric shade is not shown in Fig. 7. Moreover, since photometric shade does not resemble a vegetation endmember, best estimates of tree shade coverage in the simulated image were obtained by not combining the image tree endmember fractions

with fractions of photometric shade. Table II specifies the treatment of shade fractions for traditional unmixing methods. Table IV shows a comparison between the real errors attained by estimating GT fractions with the mean bundle images (both raw and perturbed), with the derived endmember fractions, with image endmember fractions, and with the GT mean fractions. The large errors in the image endmember fractions are partially due to the fact that in the GT fraction images, the maximum pixel fraction for trees is 0.74 and for senescent grass 0.52. That is, there are not sufficiently pure pixels of these ground components to be used as endmembers. The errors using derived endmembers are due to endmember variability, which forces the derived endmembers not to coincide with endmembers used in the simulation bundle. The mean bundle and the mean GT fractions have very similar errors. Moreover, our procedure to perturb the mean bundle fractions so that they sum to 1 did not introduce new error into fraction estimates (Table IV). Thus, the mean bundle and perturbed mean bundle methods perform as well as a complete collection of endmembers in the field (which would not be reasonable to accomplish). Moreover, with Dikin’s affine


1091

Fig. 7. Reflectance spectra of pixels with highest fractions of tree, senescent grass, and soil. These spectra were used as image endmembers to unmix the AVIRIS simulated scene.

Fig. 8. Averages of all tree, senescent grass, and shade endmembers created by the canopy RT model and combined to form the pixel reflectances in the simulated image. These spectra have been reduced to the first four principal components of the simulated data. TABLE IV AVERAGE ERROR IN FRACTIONS UNDER FIVE DIFFERENT COMPUTATIONS OR UNMIXINGS. FRACTION VALUES USED IN THE SIMULATION ARE CONSIDERED THE GROUND TRUTH. STANDARD DEVIATIONS ARE ALSO SHOWN

the GT mean tree image is somewhat lighter, implying larger tree fractions. The image endmembers had fraction images in which the spatial distribution of soil and senescent grass appear confused, probably because of the extreme variability in the two senescent grass bundles (Fig. 5). The derived endmember tree fractions underestimate the true values. The mean bundle fraction images are not shown in Fig. 9 because they are very similar to the perturbed mean bundle images. IV. CONCLUSIONS

linear programming algorithm, estimates of the error bounds accompany the fraction estimates. The images in Fig. 9 exhibit agreement and disparities among the GT mean, perturbed mean bundle, derived, and image endmember cover fractions. The GT mean and perturbed mean bundle images are similar for senescent grass and soil, but

Variation in the structural and biochemical attributes of vegetation is sufficient to cause large errors in fractional estimates computed from spectral mixture analysis. In this study, the maximum error in the tree GT mean fraction image was 0.37. Endmember bundles provide a conceptual framework for understanding this variability and for quantifying the uncertainties it introduces into fractional estimates. Averaging of the products of bundle unmixing (i.e, maximum and minimum fraction images) provides another method for estimating endmember fractions that can be significantly more accurate than methods based

1092


Fig. 9. Fraction images for image endmembers, derived endmembers, perturbed mean bundles, and mean GT endmembers.

on image endmembers, derived endmembers, or an under-representative set of field endmembers. These averages can be used to compute a new set of fractions estimating the GT fractions with approximately the same accuracy but having the additional, often-desired property of summing to 1. Finally, despite the presence of nonlinear mixing and noise, the bundle unmixing procedure found good estimates of cover fractions and their errors. These estimates could supply the accuracy of cover information needed for both the structural analyses of land use/cover change and for initializing and constraining ecosystem models. APPENDIX SIMULATED ANNEALING ALGORITHM FOR AUTOMATICALLY DERIVING ENDMEMBERS We describe a method for constructing a simplex from a partition of the facets5 of the convex hull of a data cloud in the following. The method of construction uniquely determines

d

5Let CV be the convex hull of a finite set of data points in -dimensional space. If the points in CV do not lie on a ( 1)-hyperplane, then they are the vertices of a set of ( 1)-simplexes completely covering the boundary of the convex closure of the data points and having the property that if any two simplexes in the set intersect, then they intersect along a common -face ( 1)-simplex in the set a facet of the convex hull. 1). We will call each (

d0

d0

d0

d0

k

k