A method for manual endmember selection and spectral unmixing

ELSEVIER

A Method for Manual Endmember Selection and Spectral Unmixing Ann Bateson* and Brian Curtiss* T h e number of spectraUg unique signatures needed to reproduce the statistically significant variance observed in multispectral and h!lperspectral datasets can be estimated from the eigenvalues of a principal component analysis (PCA) of the data. In this article, we describe a multidimensional T:isualization method for interactively searching for a feasible set of spectral signatures in the space of the PCA eigenvectors that account for most of the variance. These spectral signatures, referred to as endmembers, are input to a linear mixture model which can be in~erted to compute endmember abundances for each data spectrum. The visualization method discussed in this article and referred to as the manual endmember selection method (MESM) is based on lnselberg's (1985) parallel coordinate representation of multidimensional spaces. It is' novel in the .field of multidimensional visualization in that it includes not onl!t a passive representation of higher-dimensional data hut also the capabilit~f to interact with and move geometrical objects in more than 3 dimensions. The spectral shape of endmembers selected with the MESM ma!t be influenced b{t processes' such as multiple scatterintzs btj surface materials and other factors" such as illumination geometr!t that affect the signal received b!! the sensor. These processes and factors' may produce sitznificant errors' in computed endmember abundances, if not accounted for in the endmember reflectance (Roberts, 1991). The MESM is one method for obtaining endmembers that account fi~r all factor,s' and processes significantly aJfectin{3 the spectral data.

(:enter fi,r the Study of Earth from Space / CIRES, Universitx of ( :olorado, Boulder Address correspondence to C. Ann Bateson, CIRES, Unix. of Colorado, (3ampus Box 216, Boulder, CO 80309-0126. Beceived 21 Nocember 1994; recised 20 July 1995. REMOTE SENS. ENVIRON. 55:229-243 (1996) ©Elsevier Science Inc., 1996 65,5 Avenue of the Americas, New York, NY 10010

INTRODUCTION

Multispectral imagery, such as a Thematic Mapper (TM) dataset, is often viewed as a collection of hand images. Mternatively, multispectral and h~perspectral imagery can be seen as a single image with a spectrum of reflectance values associated with each image pixel. These spectra may be recognizable as the signatures of ground materials such as green vegetation, soil, or rock, provided that the material fills or almost fills the location in tile scene corresponding to the pixel. Spectra of mixtures can be analyzed with linear spectral mixture analysis (LSMA), which models each spectrum in a spectral dataset (not necessarily an image) as a linear combination of a finite number of spectrally distinct signatures, with coefficients or fractional abundances between 0 and 1 and summing to one (Adams et al., 1985; Smith et al., 1990; Roberts, 1991). These signatures are referred to as endmembers. Applied to narrow band hyperspectral data, LSMA may be preferred to band ratios such as the simple ratio (SR = R E D / N I R ) and the normalized difference vegetation index [NDVI = (NIR - RED) / (NIR + RED)] for tracking spectrally defined materials (i.e., endraembets in LSMA) since it uses all the bands in the data. LSMA is particularly nsefid when applied to an image. In this case, it produces a suite of ahundanee images, one for each endmember in the model. Like band ratio images, each abundance image shows the spatial distribution of the spectrally defined material. In addition, the spectral endmember may he identified with a specific ground material, such as green vegetation, nonphotosynthetic vegetation, or soil, by spectral f}eatures (e.g.:. the chlorophyll absorption feature at 680 nm, the water absorption fbatures at 1000 nm and 1200 nm and the NIR plateau in green vegetation) that characterize the reflectance properties of the ground materi~d at the spectral resolution of the instrument. That is, in LSMA materials can be identified as well as spatially discriminated. 0034-4257 / 96 / $15.00 SSDI 0034-4257(95)00177-:3

230

Bateson and Curtiss

Inversion of the LSMA model to compute abundances is easily performed (Boardman, 1990) once the endmembers have been supplied. In this article, we present a method that allows the researcher to construct from the spectral dataset itself the endmembers to be used in the model. The number and characteristics of these endmembers in a spectral dataset are determined not only by the number and characteristics of the spectrally unique materials on the surface but also by illumination geometry and processes (e.g., multiple scatterings by surface materials and atmospheric scattering and absorption) affecting the signal received by the sensor. A common approach for determining the number of endmembers (e.g., Smith et al., 1985) is to view the spectral data as points in the space of the complete set of spectral bands and to find the PCA eigenvectors or directions accounting for most of the variance in the data. The remaining variance is equal to instrumental error. Since sampling a continuous reflectance spectrum in many narrow, contiguous spectral bands results in a high covariance between the bands, the number of useful eigenvectors will be far less than the number of bands in hyperspectral datasets. If there are N such eigenvectors, then there are N+ 1 endmembers in the space they span, which will be referred to as the mixing space. The mean-corrected I data spectra can be orthogonally projected onto a scatter of points in the mixing space, and unmixing can be performed in that space, once the N+ 1 endmembers have been selected. The manual endmember selection method (MESM) discussed in this article is a multidimensional visualization technique for interactively exploring the mixing space in search of spectra to designate as endmembers. It relies on user intervention to judge the resemblance of spectra in the mixing space to known endmembers. In particular, selected endmembers should have nonnegative reflectance. Moreover, the endmembers should contain the data within the simplex 2 they span in order for fractional abundances to lie between 0 and 1 and sum to 1. The MESM has been implemented in a computer display, the manual endmember selection tool, which will be used in this paper to illustrate the concepts of the visualization. The MESM is novel in the field of multidimensional visualization in that it includes not only a passive view of higher dimensional data but also the capability for the user to manipulate objects in higher dimensions. Parallel coordinates, developed by Inselberg (1985), have provided a compact and structurally rich system in which to represent data in more than 3 dimensions.

~Tbat is, the average spectrum of the dataset is computed and subtracted from each data spectrum. 2That is, a triangle in the plane, a tetrahedron in 3-space and in N-space a polyhedronwith N+ 1 vertices and an edge tbr each vertex pair.

The literature includes applications of parallel coordinates to air traffic control (Inselberg and Dimsdale, 1994; Inselberg et al., 1991), robotics (Cohan and Yang, 1986), chemistry (Luke, 1993), aerospace engineering (HeUy, 1987), and statistics (Wegman, 1990). In previous applications of LSMA, selection of endmembers has been: 1. From the spectra corresponding to pixels within an image, 2. From libraries of endmember spectra (Boardman, 1990) either measured in the laboratory or in the field, 3. With a two-step method that models the spectral data with endmembers selected from the image, which are in turn modeled as mixtures of library endmembers (Smith et al., 1990; Roberts et al., 1993), 4. With automatic methods of constructing endmembers from higher order PCA eigenvectors (Full et al., 1982; Erlich and FulP 1987; Boardman, 1993). Choosing spectra fi-om the image as endmembers is perhaps the easiest method of endmember selection. However, the image may not contain sufficiently pure examples of a material (e.g., dead vegetation in an image collected during the growing season) to allow its identification from the spectra of individual pixels. With the MESM, however, endmembers need not coincide with image spectra. If a spectral endmember has snfficient variation in the image to influence the eigenvectors, the researcher can construct that endmember with the MESM. The MESM, which does not require a library, is an easier method to use than library based methods because of the difficulty in collecting a spectral library that adequately accounts fi)r all processes and factors influencing the data spectra. For example, the library mav need to include several spectra corresponding to the same green vegetation type but over different backgrounds, since multiple scatterings between leaves and a bright soil background increase the NIR reflectance of the leaves (Roberts, 1991; Borei and Gerstl, 1994). The effect of this multiple scattering on refectance varies with the mlmber of leaf layers, soil reflectance aim leaf transmittance (Borel and Gerstl, 1994; Roberts, 1991). Also, the library may need to include the effects of biological processes such as senescence and tluorescence in response to stress, since these processes can affect the chlorophyll reflectance spectrum and the spectral shape of vegetation endmembers [see Gamon et al., (1990) tbr the effects of fluorescence on vegetation +In Full et al. (1982) and Erlich and Full (1987) eudmernhcrs are not spectral, ttowever, their methods could be applied to spectral endmemhers.

Endmember Selection and Spectral Unmi~'ing 2 3 1

reflectance]. Hence, access to a very large library or extensive field work is needed for a library-based method for choosing endmembers. Moreover, when different methods of converting radiance into reflectance and removing atmospheric effects are applied to data and library spectra, often the data reflectance spectra are not accurately aligned with library reflectance speetra (Kruse, 1992; 1993). On the other hand, this problem is eliminated with the MESM, since endmembers are selected from the data space itself. Conversion techniques and atmospheric removal methods are adequate for the MESM as long as they result in spectral endmembers recognizable hy the researcher. Automatic methods of constructing endmembers fiom the eigeuvectors have the advantage of eliminating human interaction time and of producing the same results each time applied. However, MESM may be preferred over an automatic method since it allows the researcher to control which spectra, from the infinite possibilities in the mixing space, will be used to analyze the data. For example, with the MESM, the researcher can often a d j u s t the depth of the chlorophyll absorption |~ature around {~tS0 nm and the height of the NIR plateau in a vegetation eudmember or the albedo of a soil endmember. Of course, automatic methods could be combined with the MESM to initialize the location ofendmembers which the MESM could then adjust. We stress that since the mixing space can contain infinitely many spectra of varying albedo resembling a given endmember (e.g., \egetation or soil), the abtmdances computed in LSMA are relative abundances. Regardless of the method used to select endmembers, calibration to gronnd measuremeuts is tweded to convert to absolute abundances. In this article, we demonstrate the MESM with 81 field spectra of alpine tundra at the Niwot Ridge Long Term Ecological Research Site in Colorado. The spectra were collected with a Personal Spectrometer II (Analytical Spectral Devices, Inc.) with a range of 350-1050 mn at a 1.4 um resolution, Each spectrum is an average of 20 readings of approximately the same ground location, with the range restricted to 350 bands between 400 nm and 900 nm and with the optic head at a height of approximately a meter from the ground. At that height, the instantaneous field of view (IFOV) of the inst,umen! has a diameter of roughly 0.6 m. Software internal to /he spectrometer converted the raw sensor response to reflectance by correcting for dark current and dividing the result by the sensor response to a Halon reflectance standard, also corrected for dark current. The S1 spectra range over six community types4: fellfield, (h~ meadow, moist meadow, wet meadow, snow bed, and shrubs. The plots corresponding to the spectra consist of various mixtures of vegetation, litter,

~%e~' Max and \\'ebber (19~2) lot a descripti(m of tbese types.

moss, rock, soil, and lichen. Ground cover was approximated with a point classification (n = 27) of each plot and the vegetation cover averaged 62% across the 81 plots (sd = 23.7%). Two plots in the moist meadows bad close to 100% green cover and one barren plot had 100% rock cover. Litter, soil, and moss/lichen never had more than 50% coverage. To demonstrate that the MESM does not require pure spectra in the data for endmember selection, we restricted the spectral data used fi)r the selection to the 71 spectra from the Niwot Ridge dataset that represent less than 90% cower of any endmember. In addition to the field spectra, 19 field endmembers were recorded with the spectrometer held at a height that excluded other materials from the sensor's IFOV. This set of endmembers includes two soil, three rock, two lichen, three moss, three dead grass, and five green vegetation (caltha and willow) endmembers and a shaded willow endmember. In the next section of this article, the mathematical steps |or linearly nnmixing spectra with a given set of endmembers are presented in order to explain details of LSMA that influenced the design of the MESM. In the third section, endmember selection in the simple ease of 2 or 3 dimensions is described. The MESM uses parallel coordinates to represent mixing spaces with dimension exceeding 3. The fourth section has a brief tutorial on parallel coordinates and a description of their application to endmember selection. The fifth section compares the endmembers derived from the fiekl spectra with the MESM to the 19 field endmembers, The sixth section concludes with a discussion of our experience with MESM. UNMIXING P R O C E D U R E

After endmember selection, unmixing of the spectral data to compute abundances is straightibrward. If the dimensionality of the data is known to be N, then, the scene has N+ 1 endmembers, Et . . . . . E~.~, which determine an N-dimensional hyperplane. 5 The endmembers and spectral data are offset by E~ so that the hyperplane passes through the origin. The resulting space is the mixing space. An N-dimensional orthonorreal basis ~ is determined for this space, and the offset data spectra are orthogonally projected into this space. The resulting spectra and offset endmembers are expressed in terms of the N-dimensional basis. For each such offset and projected spectrum p = (2.~. . . . . )t~) in the data, solution of the matrix equation below yields the abundances (ct . . . . . c~+ ~) of the N+ 1 endmembers:

~An N-dimensional hyperplane is the set of all points in (N+ 1)dimensional space satisfying a linear equation with N + 1 variab}cs. ~The vectors in an orthonormal basis are orthogonal to each other and have trait length.

232 Bateson and Curtiss

I0

-1 C12 ... gIN+I/ • . .

~

CN

1

1

1

1 ]

C,~,+

1

The vector (e, . . . . . e~.0 is the offset eigenvector E, expressed in terms of the new basis. Note that the last row of l's guarantees that the abundances sum to 1. A spectrum p viewed as a point in N-dimensional space lies within the N-simplex with E~. . . . . E, +, as vertices if and only if 0 < c~ < 1. The endmember selection tool aids the researcher in finding endmembers as vertices of an N-simplex containing the spectral data.

LOWER-DIMENSIONAL E N D M E M B E R SELECTION

The first eight PCA eigenvectors of the Niwot Ridge field spectra are shown in Figure 1. In this article, we will use the MESM to construct six endmembers from the first five eigenvectors in this set. Corresponding to each eigenvector, the PCA yields an eigenvalue which equals the variance in the data explained by that eigenvector. The first four eigenvectors in Figure 1 account for 99.8% of the total variance. The unexplained variance is close to the total variance attributable to the instrument, computed for a signal to noise ratio of 2000: 1 and a dataset of 71 spectra, which are each an average of 20 spectral readings with 350 bands. However, since noise can occur in the higher order eigenvectors, there can be some signal in the remaining eigenvectors. Because the residuals between some field spectra and their projections into the space of the first four eigenvectors did not appear to be due to noise (see Fig. 2), we decided to select endmembers in the space of the first five eigenvectors. Also, during the endmember selection process, we obtained endmembers with unusual features that could be eliminated by selecting endmembers in a 5-dimensional mixing space (see the fifth section on Results of the Selection)• When the noise structure is known or can be estimated from a hyperspectral image, the MNF transform (Green et al., 1988) can be used to determine the dimension of the mixing space. The forward MNF transform is used to equalize the noise variance to 1 in all directions, whence eigenvectors with variance equal to 1 are all noise and can be discarded. The higher-order images can be viewed in sequence to determine if they have useful spatial information. The dimension of the dataset equals the number of useful eigenimages. Endmember selection with the MESM is currently done incrementally. That is, the MESM starts with the selection of endmembers in the plane and proceeds to higher dimensions, adding a new endmember with each dimension. For the construction of three and four endmembers, the MESM uses a simple projection method rather than parallel coordinates. Figures 3 and 4 show

the displays of the manual endmember selection tool developed for selecting three and four endmembers. In Figure 3, the E n d m e m b e r Selection Screen displays the mean-corrected spectra projected into the plane of the first two eigenvectors. As the cursor is moved over the coordinate plane, the Current Spectrum window displays the spectrum of the point underneath. Since this spectrum is formally constructed as a linear combination of the eigenvectors, it can have negative values. In looking tor endmembers, the user looks for points whose spectra resemble known signatures and hence have positive reflectances between 0 and 1. In Figure 3, three endmembers have been selected. If the dataset were 2-dimensional, the endmembers would have been selected so that all projected spectra lie within the triangle they span. When a fourth endmember is selected, the plane of Figure 3 becomes 3-dimensional space, and the fourth, automatically generated endmember is the apex of a tetrahedron whose base is the triangle constructed in the plane. The endmember selection screen in Figure 4 shows the data projected from an endmember into the plane containing the face of the tetrahedron opposite that endmember. The Rotate Tetrahedron window shows the currently constructed tetrahedron. The planes of the other faces of the tetrahedron can also be viewed in the Endmember Selection Screen. The user can continue to travel through these planes looking for spectra that are pure endmembers and relocating the vertices to enclose the projected points, if desired. Figure ,5 shows the four endmembers selected in 3-dimensional space. HIGHER-DIMEN SIONAL E N D M E M B E R SELECTION

Generalization of the projection method to arbitrary dimensions entails projecting the data onto triangular faces of the constructed N-simplex and manipulating the vertices of these faces. However, because of the following shortcomings of the projection method, a method based on parallel coordinates was developed. 1. The coordinate system becomes confusing when the data are projected onto triangular faces. That is, it is difficult to change a spectrum only in the direction of a specified eigenvector. 2. It is difficult to reach all points in higherdimensional space by relocating the vertices of the triangular faces. 3. The number of triangular faces grows on the order of N a, where N is the dimension of the space. A Brief Overview of Parallel Coordinates The polygonal line in Figure 6a is the parallel coordinate representation of the point (0, 4, 2, 5, 1) in 5-dimensional

Endmember Selection and Spectral ['nmixin~z

E|cjenvector 1 0.12

.........

, .........

, .........

, .........

, ...................

nn~]

........

Eigenvolue=

0.10

o.~

~

o-°4

, .........

, ........

7.60000

.

Eigenvolue=

O.Oe

~

L';gm~tctor 2 , ......... , .........

, .........

233

t .-0.oe

0.00 -0,,02 3OO

~

_

. . . . . . . . 4OO

~

-0.08 , ................... ~X)

, . . . . . . . . . . 7OO BOO

....

-0.1o SmO

i .........

J .....................

400

3o0

~.,(..,tn)

~Q

[igenvector 0.15

...................

0.10

3 . .........

, .........

n

700

I[;~lenvector , .........

, .........

0415

~

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

T

, .......

600

.

.

.

.

.

.

.

.

.

,

.

.

BOO

gco

4 .

.

.

.

.

.

.

]

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

o.~0 E i g e

0.0~

.

0.00

-0.05

-0.1~,

-o.lo E g e

-0.10

-0.1~

......

i ......... 400

.'TOO

~,, 500

.......................... eO0 ~,,~u,(um) [;cjenvector

O.

10

.

.

.

.

.

.

.

'

. . . . . . . . .

'

. . . . . . . . .

'

700

, ........ 800

-0.1"

......................................................... 400 ~Q e00 700

~0

5

. . . . . . . . .

8110

w,,..*,,,~.m) '

. . . . . . . . .

0.20L .........

~ . . . . . . . .

' .........

' .....

F'kjec~ector 8 ' ' ' .........

' .........

' ........

0.15~ OAO~-

0.~ -0,10

-0.06

Eige -0.15 -0.20

-0.10

300

400

~

8(10

700

BOO

E[genvalue=

-0.15

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

gO0

0.0120000

. . . . . . . . . . . . . . . . . . . . . . . . . . .

3OO

400

500

E;genvector 7 0.20

........

, .........

, .........

' .........

i .........

, ........

7100

~

ImO

o.,~, ...................... E:!I~'?.'~ .e ...................... 0.,Ol Eige

0,1~ 0.10

0.1~1

9 ~

(too W,,,,,W,,~(,~]

It~mkmgth(~q.n)

0.~1 0.~

2 0.0o6ooooo

. . . . . . . . . . . . . . . . . . . .

300

400

500

600 t l l ~ t

700

-0.1S

800

h(um)

.........

300

, ...................

400

SO0

, ............................

~

7100

800

"m,,W,~,(,,',,)

l,'i~ure 1. The first eight ei~envectors from a principal c(mq)onent analysis of tht, Niwot Rid¢e field spec'tra.

space. The five equally spaced ~ertical lines are the parallel axes of the coordinate system, and the polygonal line is the plotted point Note that spectra have a very natural representation in parallel coordinates with a parallel axis for each wavelength. A PCA of the spectral data may be used to reduce the dimension and, hence, the number of parallel axes While points are represented 1)v polygonal lines

ill parallel coordinates, lines are duallx represented ~s points. Figure 61) shows points in parallel coordinates lying on the line !t = m*x + b. The point 1 is said to represent the line ¢!= m*x + 5 in parallel coordinates. Observe that if l~ and 12represent parallel lines, then they lie on the same vertical line in parallel coordinate representation. Figure 7a shows the location of points representing lines in parallel coordi-

234

Bateson and Curtiss

First 0.3

. . . . . . . .

I . . . . . . . . .

i

4

Eigenvectors

. . . . . . . . .

i

. . . . . . . . .

i . . . . . . . . . . . . . . . . . . .

0.2 L)

~- 0.1

0.0 300

400

500

600 700 Wove'ength

800

900

First 5 [igenvectors 0.3 . . . . . . . . . ' . . . . . . . . . . . . . . . . . . . . . . . . . . .

®

0.2

C

o~ 0.1

0.0 300

400

500

600 700 Wave!ength

800

900

Figure 2. Comparison of a field spectrum from a sparsely vegetated plot with its image in the space of the first four eigenvectors and the space of the first five eigenvectors. The field spectrum is plotted with a dashed line and the image spectrum with a solid line.

nares as a function of slope. The extension of parallel coordinate representation of lines in the plane to higher dimensions is illustrated in Figure 7b, which is a parallel coordinate plot of collinear points in 5-dimensional space. The parallel coordinate treatment of planes is based on the simple observation that two intersecting lines determine a plane. In the Cartesian representation in Figure 8a, the lines are l~ and 12, which are the intersection of the plane with the xy and yz coordinate planes, respectively. Given point P on the plane containing It and l~, the point W which lies on the intersection of lz and the line through P parallel to line ll has the same z-value as F. The above observations are used in Figure 8b to illustrate how the z-value of the point P can be determined graphically in parallel coordinates from its x- and

y-values. The polygonal line labeled with x', g', and z' represents the point P whose z-value we are computing with the aid of the point W (0, y", z'). Note that the parallel coordinate representations of P and W cross the vertical line through the point l~ (representing in parallel coordinates the line l~ of Fig. 8a) at the same point, since W lies on the line through P parallel to l~, Also, the parallel coordinate representations of A and W intersect at the point 12 (representing line l~ in Fig. 8a), since they both lie on the line l~. These two observations, together with the fact that W lies in the yz-plane, determines W and hence z'. In higher-dimensional endmember selection, the lace of an N-simplex can be translated away f~om the vertex (endmember) opposite that face by adjusting A. The face can be rotated by adjusting the location of the vertical lines through l~ and l_o.

Application of Parallel Coordinates to Endmember Selection For endmember selection ill a mixing space with more than 3 dimensions, the MESM has three basic parallel coordinate representations which allow the user to: 1) relocate individual endmembers, 2) determine if spectra lie outside of the current endmember simplex, and 3) translate and rotate faces of the simplex to encompass spectra that lie outside the simplex. The manual endmember selection tool implements these three representations in the Data, Cluster, and Face windows respectively. These windows, showing the Niwot Ridge field data in the mixing space of the first four eigenvectors, are seen in Figures 9 and 10 and discussed in detail below. 1. Data Window: The data window in Figure 9a has a parallel axis for each eigenvector in the mixing space. The black polygonal lines are plots of the data spectra orthogonally projected into the space of the first four eigenvectors, and expressed in terms of these eigenvectors, which torm an orthonormal basis for the space. The colored lines are plots of the current endmembers. The first four endmembers are those selected with the 3-dimensional display, and the last endmember was automatically selected by the tool. Coordinate values for the endmembers can be adjusted in the data window. 2. Cluster Window: The cluster window in Figure 9b has a parallel axis for each endmember. The endmember abundances for each spectrum are plotted on these axes. The user can look for and select clusters in this window. For example, the user could select the cluster of all spectra which have a negative abundance of the first endmember E~. These spectra lie outside of the simplex on the wrong side of the face opposite El. Fig-

Endmember Selection and Spectral Unmixing 2 3 5

~.60

0,50 0,40

0~0

/

0.20

J

0.10 O.oo 300

500

4430

6(El 700 wc2w~ngth (~m)

l i t

~ l l f r ,

800

~30

I

.m

Ey

[3

*:

~.:::~

_4_, ,

-2

0

2

a

6

Figure 3. The endmember s('-

lection tool with three selected endmembers.

ure 9c shows the effect of selecting all spectra on the wrong side ot" the face opposite E4. 3. Face Window: While the data window allows one to adjust the location of individual endmemhers, the face window allows the user to rotate and translate whole (N-1)-dimensional faces. Figure lOa shows the t:ace window after the selection of all spectra on the wrong side of the face opposite E~. These spectra are displayed in the t:,tce window in the following way. K the equation of the hyperplane is (!I'XI

+C2"Xi+C2"X3+C4"X4

~

CO,

then the face window includes the vertical lines through the points It and 12 representing lines It and 12, as in Figure 8a, to represent the plane C2"X + C3" y + C4" 7, = 8 0 .

In Figure 10a, the green and pnrple vertical lines pass

through the points l~ and/,2, respectively. Each spectrum in the cluster is mapped onto a point (x,y,z) in the space of the face window by the transformation x = (cl "Xl + c~'x2)/c~, y = X3, ,Z = X~ .

We will assume the data is translated so that point A, where the lines represented by 1~ and 12 intersect, coincides with an endmember. In Figure 10a, A is the magenta polygonal line. Each transformed spectrum (a,b,c) in the cluster is displayed in the Face window along with the helping point that computes the value c' that would place (a,b,c 3 on the plane. In Figure 10a, data spectra are represented with black polygonal lines and the helping points with blue polygonal lines. Since all the spectra in the cluster lie on the wrong side of the face F, the value of c is

236

Bateson and Curtiss

m

Figure 4. The endmember se-

lection tool after a fourth endmember has been selected.

either too small for all the selected spectra or too large. By adjusting the y-coordinate of A, the user translates the hyperplane away from E4. By observing the changing relationship between the helping points and the selected spectra in the cluster, the user can see when the hyperplane has been translated beyond the spectra in the cluster so that they now lie on the right side of F--i.e., they now have positive abundances of endmember E4. Note that this translation changes all endmembers but E4. As an alternative, the user can rotate the hyperplane by adjusting the locations of the vertical lines through Ii and 12. Figure 10b shows the relative position of the spectra and helping points after the user has moved the purple vertical line through 12 to the right. In summary, since the face window manipulates whole (N-1)-dimensional faces instead of the vertices of 2-dimensional

laces, the complexity of the display is reduced from N :~ to N.

RESULTS O F THE SELECTION

The final endmembers selected in 4-dimensional space are in Figure 11 and visually identified as vegetation, shade, dead vegetation, rock, and soil. W e continued endmember selection in 5-dimensional space to improve endmember fits to spectra containing large rock or soil abundances and to eliminate the unusual protuberances in the rock and soil endmembers centered around 720 nm. As seen in the cluster window (Fig. 12 ), after the six endmembers appearing in Fig. 13 were selected, we allowed some small negative abundances, which we considered to be 0.

Endmember Selection and Spectral Unmixing

23 7

0.6

{3,4 2

o

fr

0,2

0.0

Figure 5. The Saved Endmemher window with four endmembers selected from the Niwot Ridge field spectra the 3-dimensional endmember selection display.

400

500

Figure 6. Illustrations of basic concepts in parallel coordinates representations: a) a point in 5-space; b) a line in the plane, 5,

3

xl

x4

(a)

b/(1 -m)

Ix

1/(l-m)[

Y

1 (b) l

600

700

way+re+ore (urn)

I

x5

800

When multiples of the fifth eigenvector were added to the rock and soil e n d m e m b e r s in Figure 11, the protuberances centered around 720 nm were smoothed out (see Fig. 13). Also, a new vegetation e n d m e m b e r was selected that has an uncharacteristic peak around 720 nm. We speculate that green vegetation fluorescence, which increases reflectance around 720-740 nm (Schmuek, 1990) was incorrectly expressed in the soil and rock e n d m e m b e r s in the space of the first four eigenvectors. However, in the 5-dimensional mixing space, another vegetation e n d m e m b e r was added that correctly accounts for fluorescence in its reflectance spectrum. The shade e n d m e m b e r in Figure 13 is interpreted as a functional e n d m e m b e r that does not represent the spectral signature of any material in the scene but, instead, accounts for the effects of shade. Note that its spectrum rises in the blue wavelengths, as we would expect, because of increased scattering in those wavelengths. The curve-fitting procedure of Gao-Goetz (1994), which compensates for albedo differences and differences in feature depths, was used to compare all endmembers in Figure 12 except for shade to the 19 field endmembers. For the comparison, the e n d m e m b e r

238 Bateson and Curtiss

(a)

/

m=

0
m<0

m>l

12

x

Y

Y

x

|

l

w

(0~',,z;)

A= (0, y, O)

(b)

- - p Y

\>

. . . .

W

=...A 12

I~s

I

/ / " y'

Xr /

% y

xl x2 x3 x4 x5 Figure 7. Additional information on parallel coordinates rep-

resentation of lines: a) location of a point representing a line as a function of the slope m; b) four polyglonal lines representing four points in parallel coordinates that lie on a line in 5-space.

spectra are divided into intervals and the curve-fitting procedure is applied to each interval. This procedure first converts reflectance into absorbances by taking the negative of the log of both the field endmember and the derived endmember to which it is being fitted. A line is removed from each endmember to remove albedo effects and a scaling factor is determined by a linear regression of the derived endmember absorbance against that of the field endmember. The best fits are shown in Figure 14. The statistics for each fit includes the scaling factor, which should be positive for good fits, and the correlation coefficient from the regression between the derived and field endmembers. The fitted spectrum is obtained from the scaled field absorbances by replacing the background line and converting to reflectance. We also tried selecting six endmembers from the library of field endmembers. In particular, we used a program that cycles through all six element subsets of the 19 field endmembers and computes for each subset the number of field spectra that have positive abundances for all endmembers in the set. The maximum

% z

Figure 8. Depiction of the concepts behind parallel coordinates representation of planes: a) a plane represented in a Cartesian coordinate system; b) illustration of how to construct the z-coordinate of a point P on the plane a) from its x and y coordinates using parallel coordinates representations of points and lines.

number of field spectra with positive abundances over all six element subsets was only 34 out of 81 with negative abundances as extreme as -0.95. This performance illustrates the difficulty in collecting a spectral library that accounts for all processes (e.g., fluorescence) 'affecting the signal received by the sensor. Finally, we computed Pearson's correlation coefficient to measure the strength of the linear relationship between the sum of the fractional abundances of the two vegetation endmembers and the ground approximation of vegetation abundance mentioned in the introduction. Over all 81 plots, Pearson's correlation coefficient is 0.71 (p < 0.0001). EXPERIENCE WITH THE METHOD

The MESM has been used to derive endmember spectra from AVIRIS (Airborne Visible and Infrared Imaging Spectrometer) data and field data collected using a Personal Spectrometer II. In Wessman et al. (1995), litter and vegetation abundances from a spectral unmixing of AVIRIS data were significantly correlated with ground measures of litter and vegetation at the Konza

Endmember Selection and Spectral Unmixing 239

C

(a) Figure 9. The Data and Cluster window of the higher dimensional display: a) the Data window on entry to the display; b) the Cluster window on entry; c) the Cluster window after selecting to view all spectra with negative abundances of the fourth endmember (the green axis in Fig. 11).

i

lal

i i ~!

!',

'

1

1l

~

w

'

~

Figure 10. The Face window of the higher-dimensional display: a) the Face window showing all spectra on the wrong side of the hyperplane opposite the fourth endmember; b) the Face window after adjusting the purple vertical line in (a) to rotate the face opposite the fourth endmember so that the spectra now lie on the same side of that face as the fourth endmember. Parallel axes represent eigenvector directions selected by the researcher with the buttons underneath the plot.

..... ir

"I

Ill

''

i~i~~ ~ ~i~ ~i~

'~ii~i~I ~

!~!i}i¸/:

i i!~'~' i!|~iill i~ i

/ilil iiilli

/

ii! (b)

rl

71ii ¸¸¸¸¸¸

?

I

ii

IOt/)

240 Bateson and Curtiss

Vegetation

Live

10 I ~ 0.4

.

.

200 400 600 8001000

Wavelength

°j i

Wovelength

0.4

0.2 0.0

08

~

06

~

-

0.4

Wovelength

Dead Veqe:a:~on !.0[

0.8

0.8 i

o6 i

0.6!

02 ~ 0el 0 2 - ..........

021 0.0 i

~

200 400 600 80@ C30 Waveiength

o2L .............. 200400 6C0 BOS~00,3 Wcve[e~gt ~

Soil

Live v e g e : c : o r 2

Rock 1.0[

08 06

• u

0.8 t

f r

06 i

~

:

041

OA

0.2 0.0 200 400 600 8001000

Shade

o4i

Wavelength

Wovelength

0.6

0.6

8

200 400 600 8001000

200 400 600 8001000

0.8

O.B

1.0!

00[ O2

.~-.-~.,,,

Soil

Rock 1.o

-

-

0.2 0.0 ~ 2004006008001000

0.2 0.0

8

g o.6

g o.6

~' 0.4

oi 1 0i/

0.8

0.8

0.8

Live Vegetation

Dead Vegetation

Shade

0 02L 2904006008001D00 Wovelength

600

200 400 8001000 Wavelength

3,0 i

o2i 200 400 608800!00C Wave!ength

-0,2L . . . . . . . . . . . 200400 600 BOO]30C Wavele~g~"

Figure 11. The five endmembers selected in 4-dimensional space with the manual endmember selection tool.

Figure 13. The six endmembers selected in 5-dimensional space with the manual endmember selection tool.

Research Area outside of Manhattan, Kansas. Significant correlations were also found between the vegetation endmember abundances computed for field spectra collected from the Konza Research Area in 1989 with a Personal Spectrometer II and field measures of APAR.

The field data from the Konza were collected in a growing season following a drought year, and each pixel had at least a 30% exposed background according to C. A. Wessman (personal communication). Thus, the MESM can derive endmembers from highly mixed data.

Figure 12. The cluster window after six endmembers have been selected in 5-dimensional space.

Cluster Window Q.I!

0.6

@.@ -@.~ 0

Restrict

!

•

3:

4

@

Endmember Selection and Spectral Unmixing 241

Library EnOmember jy92wi13.mn 0.50

Library Endmember jy921it2.rnn ....... , .................r ' ~ "

0.20 ~

S

040 030

0.10

020

0.05 0.10 000 300

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

400

500 600 700 Wavelenglh (urn}

800

900

Endmemberjyg2wi13mn Fittedto Live Vegl Endmember 0.8

...........

' ......... ' ......... ' ......... ' ........ ; t 1

.... J...............................

0.00 300

400

500 600 700 Wavelength (urn}

800

go0

Endrnemberiy921i12.mn Fitted to Dead VegetationEndmember 0.08 ..................

, . . . . . . . .

,

............... - /

00,I

y

/

g ~ 0.4

0.04

. . /

0.2

0.02

2

o.o

Figure 14. Curve-fits of field endmembers to five of the six endmemhers in figure 12. Intervals of the eur~'e-fit are shown in the plots of the field endmemhers. Dashed lines are resealed library eudmembers.

..............

300

400

................. !

500 600 700 Wavelength (urn}

800

900

interval staUsacs for Fit of image willow endmember with library endmember willow O. [ 393.400, 731.9001 Scaling: 1.97671 Correlation: 0.928827 RMS: 0.0144386 1. [ 731.900, 877.6001 Scaling: 1.83063 Correlation: 0.995934 RMS: 0.00315"/62 "Iotal RMS = 0.0121698

0.00 ........ ~ 300 400

....................... ................. SO0 600 700 800 900

Waveleagth(urn)

Interval statistics for Fit of image dead vegetation cndmmnber with library ¢ndmember ~iRcr 0. [ 393.400, 87%6001 ~ g : 3.56025 Correlation: 0.850612 RMS: 0.00275555 Total RMS = 0.00275555

Library Endmember jyg2soi2.mn

Library Endmember jy92rockmn 0.250.20 .......; 0.20 015

,.t 0.15

o,o .................................................. i li 0.05 300

400

500 600 700 Wavelength (urn)

800

900

Endmemberjy92rock,mn Fitted to Rock Endmember 0.50

..................' ..................' .................

0.40

0.10 0.05 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30O 400 5OO 60O 700 Wavelength (urn)

8OO

90O

Endmemberjy92soi2.mn Fitted to Soil E n ~

°35f ................................................... i.;

g~ o.30 0.20 ~. 0.20 0.10 0.00 ................................................. 300 400 500 600 700 800 900 Wavelength (urn)

Figure 14. Continued.

interval statistics for Fit of image rock endmember with library endmember rock 0. [ 393.400, 877.600] Scaling: 2,73087 Corcelation: 0.956314 RMS: 0.0196070 Total RMS = 0.0196070

ZI--/ 30O

4O0

5O0 600 70O Wavelength (urn)

800

9O0

Interval statistics for Fit of imagcndmcmber soft with library endmember soil O. [ 393.400, 877,600] Scaling: 0.670359 Correlation: 0.975231 RMS: 0.00459183 Total RMS = 0.00459183

242 Bateson and Curtiss

Library Endmember jyg2clth.mn 0.8! ....... 0.6

0.4-

0.2-

0.0'. ...... 3O0

~ . . . , ~

................

50O 6O0 700 8OO 90O Wavekmgth (urn) Endmend~r jy92¢lth.mn Fitted to Live Veg2 Endmember 1.0

400

.........

' .........

' .........

' .........

~

'

'

~

'

"

/

0.8 0.6 0.4

0.2 O.O

300

.... .

........................................

400

500

600

Wavelength

700 (urn)

800

endmember in the mixing space must be relocated in order for the data to lie within the simplex the endmembers span. Preliminary endmember selections with the library link greatly reduced human interaction time. Currently, the user must select endmembers incrementally, starting with the selection of three endmembers in 2-dimensional space and proceeding to higher dimensions, adding a new e n d m e m b e r with each dimension. Future work should include combining the MESM with methods that automatically generate a set of endmembers in the final mixing space which the MESM could then refine. Such a combination would not only reduce human interaction time but also help to delimit the endmembers selected by different users. An important advantage of the MESM is the capability of deriving endmembers whose spectra reflect some of the scene-specific processes affecting the signal received by the sensor. With datasets including endmembers influenced by such processes (e.g., Wessman et al., 1995), the MESM was an important exploratory device for determining the potential of LSMA to track scene components across a landscape.

900

Intervalmtiaties for fit o f i m a g e e m d m ~ b e r vegetation2with library~ndmeml~rcdthL 0. [ 393.400,694.500] Scaling: 1.14043Con~lation:0.928928 RMS: 0.0145390 1. [ 694.500,877.600]Scaling:0.950098Correlation: 0.986906RMS:0.0330659 TotalRMS= 0.0233244 Figure 14. Continued.

A limitation of the method is that sufficient variation in abundances is required to determine eigenvectors that distinguish between the endmembers. For example, if, in the extreme, soil and litter occur in the same proportion relative to each other in all pixels, then the MESM will find a litter/soil e n d m e m b e r rather than two distinct endmembers. The weaknesses of the method are: 1. It requires human intervention time to judge the acceptability of the endmembers. 2. Different users could find somewhat different endmembers. A link between the endmember selection tool and spectral libraries is a recently added feature to give the user more direction. This link entails projecting library spectra into the eigenvector space of the data. Various quantities are computed to determine how well each library spectrum is modeled in that space and how effective each spectrum would be as an added endmember or as a substitute for a current endmember in the mixing space. Usually the projected image of a library

This research was supported by NASA/ EOS Grant NAGW 2662. Data were provided by the Niwot Ridge Long Term Ecological Research Project (NSF DEB 92117760) and the Mountain Research Station (BIR 9115097). We give special thanks to Carol Wessman for the use of the data and to Elizabeth Nel for her care and diligence in performing the spectral measurements.

REFERENCES Adams, J. B., Smith, M. O., and Johnson, P. E. (1985), Spectral mixture modeling: a new analysis of rock and soil types at the Viking Lander I site, J. Geophys. Res. 91:8098-8112. Boardman, J. (1990), Inversion of high spectral resolution data, Proc. SPIE Tech. Symp. Opt. Electro-Opt. Sensors 1298:222-233. Boardman, J. (1993), Automated spectral unmixing of AVIRIS data using convex geometry concepts, In Summaries of the Fourth Annual JPL Airborne Geoscience Workshop, pp. 1114. Borel, C. C., and Gerstl, A. W. (1994), Nonlinear spectral mixing models for vegetative and soil surfaces, Remote Sens. Environ. 47:403-416. Cohan, S. M., and Yang, D. C. H. (1986), Mobility analysis of planar four-bar mechanisms through the parallel coordinate system, J. Mech. Mach. 21:63-72. Ehrlich, R., and Full, W. E. (1987), Sorting out geologyunmixing mixtures, in Use and Abuse of Statistical Methods in Earth Sciences, W. B. Size, Ed.), Oxford University Press, New York, pp. 33-46. Full, W. E., Ehrlich, R., and Bezdek, J. C. (1982), Fuzzy Qmodel- a new model approach for linear unmixing, Math. Geol. 14(3):259-270. Gao, B. C., and Goetz, A. F. H. (1994), Extraction of dry

Endmember Selection and Spectral Unmixing

leaf spectral features from reflectance spectra of green vegetation, Remote Sens. Environ. 47(3):369-374. Gamon, J. A., Field, C. B., and Bilger, W. (1990), Remote sensing of the xanthophyll cycle and chlorophyll fluorescence in sunflower leaves and canopies, Oecologia 85(1): 1-27. Goetz, A. F. H., Vane, G., Soloman, J. E., and Rock, B. N. (1985), Imaging spectrometry for earth remote sensing, Science 228:1147-1153. Green, A. A., Berman, M., Switzer, P., and Craig, M. D. (1988), A transformation for ordering multispectral data in terms of image quality with implications for noise removal, IEEE Trans. Geosci. Remote Sens. 65-74. Helly, J. J. (1987), Applications of parallel coordinates to complex system design and operation, in Proc. National Computer Graphics Association 1987, Fairfax, VA, Vol. III, pp. 541-546. Kruse, F. A. (1992), Prototyping HIRIS analysis techniques using AVIRIS data, Final Report, U.S. Geological Cooperative Agreement 14-08-0001-A0746. Kruse, F. A. (1993), The effects of AVIRIS atmospheric calibration methodology on identification and quantitative mapping of surface mineralogy, Drums Mountains, Utah, in Summaries of the Fourth Airborne Geoscience Workshop, pp. 101-104. Inselberg, A. (1985), The plane with parallel coordinates, Special Issue on Computational Geometry of the Visual Computer 1:69-97. Inselberg, A., Boz, M., and Dimsdale, B. (1991), Planar conflict resolution algorithm for air traffic control and the one-shot problem, IBM PASC Teeh. Report G320-3559, Palo Alto, CA. Inselberg, A., and Dimsdale, B. (1994), Multidimensional lines If: proximity and applications, SIAM J. Appl. Math. 54(2): 578-596.

243

Luke, B. T. (1993), Displaying the structure of molecules by multidimensional plots of their torsion angles, J. Chem. Info. Comput. Sci. 33:135-142. May, D. E., and Webber, P. J. (1982), Spatial and temporal variation of the vegetation and its productivity on Niwot Ridge, Colorado, in Ecological Studies in the Colorado Alpine (J. C. Halfpenny, Ed.), University of Colorado, Institute of Arctic and Alpine Research, Boulder, pp. 35-62. Roberts, D. A. (1991), Separating spectral mixtures of vegetation and soils, Ph.D. thesis, University of Washington, June. Roberts, D. A., Smith, M. O., and Adams, J. B. (1993), Green vegetation, nonphotosynthetic vegetation, and soils in AVIRIS data, Remote Sens. Environ. 44:255-269. Schmuck, G. (1990), Applications of in vivo chlorophyll fluorescence in forest decline research, Int. J. Remote Sens. 11(7):1165-1177. Smith, M. O., Johnson, P. E., and Adams, J. B. (1985), Quantitative determination of mineral types from reflectance spectra using principal components analysis, Proc. Fifteenth Lunar Planetary Sci. Conf Part 2, J. Geophys. Res. 90: C797-C804. Smith, M. O., Ustin, S. L., Adams, J. B., and Gillespie, A. R. (1990), Vegetation in deserts: I. A regional measure of abundance from multispectral images, Remote Sens. Environ. 31:1-26. Smith, M. O., Adams, J. B., and Sabol, D. E. (1994), Spectral mixture analysis-new strategies for the analysis of multispectral data, in Imaging Spectrometry-a Tool for Environmental Observations (J. Hill and J. Megier, Eds.), Brussels, pp. 125-143. Wegman, E. (1990), Hyq3erdimensinnal data analysis using parallel coordinates, J. Am. Stat. Assoc. 85:664-675. Wessman, C. A., Bateson, C. A., and Benning, T. (1995), Detecting fire and grazing patterns in a tallgrass prairie using spectral mixtnre analysis and NDVI, Ecol. Appl., in review.