A Landsat MSS time series model and its applications in geological mapping

PHOTOGRAMMETRY 81REMOTE SENSING ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

A Landsat MSS time series model and its application in geological mapping Ding Yuan *, James R. Lucas, Donald E. Holland NASA Commercial Remote Sensing Program, Lockheed Martin Stennis Operations, John C. Stennis Space Centel; Mississippi, USA

Received 5 May 1997; accepted 24 July 1997

Abstract Although temporal factors in remotely sensed data have long been recognized by remote sensing scientists, studies on the separation of temporal and stable factors have been very limited. Applications of the time-variant characteristics of the data are also limited in the environmental monitoring and change detection fields. The merits and applications of the stable factors in the time-variant historical image data or image time series have not been fully recognized. This paper presents a time-variant model for the historical multitemporal multispectral scanner (MSS) data and a technique for processing the image time series for geological mapping. This model has been successfully applied to the surface geological mapping of a region in central eastern Nevada, USA. 0 1998 Elsevier Science B.V. All rights reserved. Keywords: geological

mapping; image time series; component

1. Introduction

Modem satellite remote sensing began in 1972 when the United States launched its first Earth Resources Technology Satellite (ERTS, later renamed Landsat). Since then, numerous satellite images of the Earth’s surface have been acquired by agencies and industries worldwide. The technology for Earth land-surface remote sensing has evolved from the earlier multispectral scanner (MSS), thematic mapper (TM), and European systeme probatoire d’observation de la Terre (SPOT) to the recently planned hyperspectral imager for the new-generation small satellite, Lewis. The new-generation data with higher spectral and spatial resolutions will soon *Corresponding author. Tel.: +l (601) 688-2905; Fax: +l (601) 688-3838; E-mail: [email protected]

separation;

remote sensing; Landsat MSS

appear on the market. This trend represents rapid progress in remote sensing technology and reflects endless searching for new data by the remote sensing community. However, rapid progress in remote sensing technology does not imply that the first-generation satellite remote sensing data, such as MSS and TM, will soon be obsolete. Rather, the first-generation remote sensing data play vital roles in environmental monitoring and global change research. Many global change research programs are inevitably resorting to the use of historical remote sensing data. For instance, the Environmental Protection Agency’s North American Landscape Characterization project was one of the most ambitious projects that involved the use of historical MSS data for documenting landscape and landcover changes for the North American continent (Lunetta et al., 1993).

0924-2716/98/$19.00 0 1998 Elsevier Science B.V. All rights reserved. PII SO924-2716(97)00027-O

40

D. Yian et al./ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

Of course, the temporal-variant characteristics, or the witness of history, do not represent the sole value of the historical remote sensing data. The first-generation remote sensing data preserved information about surface geology, soil types, and cultural features. Intensive human activities are continually destroying or modifying the surfaces from which this information was taken. Historical image data are better than the new data for mapping geology or soil in areas of great surface disturbances. For example, the U.S. Geological Survey (USGS) is currently mapping the geology of a 7.5~minute quadrangle in Clark County, Nevada. This area has undergone extensive development in the past several decades due to the booming gaming industry in Las Vegas. Historical air photos and Landsat MSS data are used extensively to recover the surface geology and lithology of the area. This paper presents a conceptual image time series model along with discussions of the characteristics, separation, and possible applications of its random, stable, and temporal components. A procedure for effectively extracting geological information from multitemporal historical MSS data for litbological mapping is discussed. This procedure will be applied to a testing area in central eastern Nevada for surface lithological mapping. 2. Review of previous work For any particular time, signals recorded in the remote sensing data (digital number or DN) are affected by the reflectance and scattering of the ground targets and the scattering and absorption of the atmosphere. Since the sensor, illumination, atmosphere, climate, phenology, and ground conditions are constantly changing, the digital numbers and the converted reflectance recorded in the remote sensing data from the same geographic location are also changing constantly. This temporal fluctuating nature of the data affects the results and accuracy of the analysis and interpretation of the remotely sensed data. Duggin (1985) summarized both systematic and random factors limiting the discrimination and quantification of remotely sensed radiance. The systematic factor can be corrected, but the random factors cannot be corrected. Duggin et al. (1985) further

showed that the temporal variations in TM data could be large enough to affect the analysis and interpretation of the data. The temporal fluctuation characteristics of multitemporal image time series have been studied by many authors from different aspects. Robinove and Chavez (1978) studied the scene-wise temporal change of albedos of more than 30 MSS scenes from central eastern Nevada and pointed out that the albedos fluctuated seasonally. However, due to the data volume and computational limitations on image time series analysis at the time, no further studies concerning the temporal characteristics of the data were reported from the authors. Yuan and Robinson (1990) and Yuan (1991) examined the seasonal effect of vegetation on lithologies in MSS time series. They showed statistically that the radiance from any given lithologies in an area are different in different seasons. Lodwick (1979) applied principal component analysis to the composite of two MSS images to estimate the ecological change between the two dates. Pickup et al. (1993) used multitemporal MSS data for delineating vegetation cover change. Due to their smaller data volume and higher frequency of acquisition, the advanced very high resolution radiometer (AVHRR) time series have been used extensively for vegetation monitoring or change detection studies. Justice et al. (1985) and Townshend and Justice (1986) studied the vegetation dynamics at the continental or global scale derived from the AVHRR time series. Malingreau and Belward (1992) studied scale factors in analyzing the AVHRR image time series in vegetation monitoring. Tappan et al. (1992) used AVHRR data to monitor the rangeland dynamic conditions in Senegal, Africa. Quarmby et al. (1993) used the normalized difference vegetation index (NDVI) time series derived from AVHRR data to establish a crop yield prediction model. However, few studies have been conducted using MSS, TM, or other high spatial, multispectral time series to solve large-scale local problems. Most current image time series studies concentrate on observing and interpreting the series. At most, analysis is performed for the numerical time series of limited parameters derived from a corresponding image time series, such as a scene- or regional-averaged NDVI, land cover type average reflectance, etc. Treatments of an image time series

D. Yuan et al. /ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

with spatial and spectral considerations are limited. For MSS, TM, SPOT, or other land observation data sets, few studies have used an image time series. The reason is obvious: there has been virtually no dependable technique or mutual image processing system that could handle the complexity of an image time series of MSS, TM, or SPGT data sets. Most image time series studies focused on the exploitation of the dynamic or the time-variant components of the time series. Few studies or applications have investigated the time-invariant factors or the stable components of the image time series, even for geological applications. Yuan (1990,199l) and Yuan and Robinson (1990) studied the MSS response of a time series from a number of surface geological units in central eastern Nevada. They concluded that responses from those units were seasonally different, implying that different criteria would be needed to interpret geology if different MSS data acquired in different seasons were used. Direct analysis of a multitemporal image time series requires that the data be originally sorted in the spectral dimension. Converting multichannel data to a single channel makes the analysis of an image time series much easier. For vegetation studies, NDVI and other vegetation indices have been used extensively to convert a multichannel AVHRR time series into a single-channel image time series (Justice et al., 1985; Townshend and Justice, 1986; Malingreau and Belward, 1992; Tappan et al., 1992). For the geological study, Yuan (1990) used albedo transformations to convert a multichannel MSS time series to a single-channel albedo time series. Given the current computational and system limitations, this transformation has proven to be effective for

41

geological mapping in non-vegetated regions. The procedure and techniques reported in this paper are a continuation of previous exploratory work in using the multitemporal MSS image time series by Yuan (1990, 1991), and Yuan et al. (1996, 1997). 3. The multitemporal

model

Conceptually an image time series consists of stable (S) or invariant, temporal (T) or variant, and random (R) or stochastic components in analogy to the conventional numerical time series (Fig. 1). These three factors affect the digital number records in the images. Intuitively, the conventional time series analysis techniques can be used for analyzing the image time series. The DN value at pixel level in the data can be expressed as a function of R, S and T components (Fig. 2): DN = F(S, T, R), where S = stable components, T = temporal components, and R = random components. In tbe simplest case, F can be assumed to be an additive function such as the one shown in Fig. 1. However, due to the spatial correlation of the pixels and the atmospheric scattering effects, the DN value for a pixel at any given time not only correlates to previous values of the same pixel, but also correlates to the previous values and the current values of the surrounding pixels. Therefore, the conventional time series may not be adequate for describing the image time series obtained by remote sensing techniques. In other words, the image time series obtained through remote sensing is not a simple sum of the individual time series at pixel level. Instead,

DN Or

Reflectance

Recorded Signal = R+S+T

fi

S = Stable Component: soil, geology, etc. R = Random Component: atmosphere, etc. T = Temporal Component: vegetation, land use, etc. r

Time Fig. 1. Schematic diagram showing the separation of recorded signal (DN) into random (R), stable (S), and temporal (T) components.

42

D. ban et al. /ISPRS Journd of Photogmmmetry & Remote Sensing 53 (1998) 39-53

Atmospheric Scattering R Components

Soil, Geology

Vegetation, Landcover change

Image time series T Components

Fig. 2. Schematic diagram showing how the RST components at image level can be used for different application purposes.

it has interactions at pixel, neighborhood, and even scene level. Fig. 2 diagrammatically shows that an image time series consists of variant and invariant components of different sources. Obviously, these variant and invariant components can be used for different purposes. For instance, the temporal component may be used for change detection and the stable component may be used for geological mapping. This paper focuses on designing and experimenting a group of procedures for separating the R-S-T components of an image time series and using the stable components for lithological mapping. The complete procedures for this research are outlined in Fig. 3. These procedures can be categorized into pre-, main-, and post-processing categories. Our focus is in the procedures in the main-processing category. In the pre-processing, a multiple image coverage of the same study area was first selected. Then, the images were registered to a digitized map, and a common area simultaneously covered by the selected image time series was identified. An appropriate multitemporal normalization was performed to enable image comparison on the same bases. Finally, an appropriate subset was selected from the study area for the surface lithological units. This subset will be used for classifiers’ training and accuracy assessment. These pre-processing procedures together produce a well-trimmed, coregistered, and normal-

ized image time series ready for analysis. The main processing procedures of the analysis consisted of the four steps given in Sections 3.1, 3.2, 3.3 and 3.4. 3.1. Albedo conversion - data reduction in spectral dimension ThesetZ = {Ii, . . .. ZK} forms a multispectral image time series for the given region. There are many ways to reduce a multichannel image time series into a single-channel image time series. In vegetation or forestry monitoring applications, vegetation indexes are often adopted. The albedo or reflectivity transformation was selected for the geological application. Let Ak = ALB(Zk) be the transformation that converts the multichannel image Zkinto a single-channel albedo image Ak. Then A = {Al, . . ., AK} is a oneband image time series. This transformation substantially reduces the amount of data used in further analysis with minimal loss of geological information. 3.2. Principal component analysis - component separation Principal component analysis has been used conventionally for data reduction. However, in this study it was used for component separation. Principal component analysis was performed on the albedo

D. Yuanet al. /ISPRS Joumal of Photogrammetly& RemoteSensing 53 (1998) 39-53

Main proceseing procedures

Proproceasing procedures

Spectral Albedo

Multitemporal

Registration

Data

Assembly

Rectification

I I

postprocessing procedures

Reduction Conversion

-

I

Classification Result

Temporal Reduction Principal Component Analysis

Separation of R, Components

S,

Geological I

T

I Production Geological

Geological Classification Using S Components

I

I

of Maps .:

I

Selection

Fig. 3. Procedures for geological mapping using multitemporal image data. The main procedures are the focal points of this paper.

time series A. The component separation ability of

the principal component analysis over multitemporal data was first observed by Ingerbritsen and Lyon (1985) for two-date composite MSS data. The application of principal component analysis to multitemporal data was tested by Yuan (1990). Generally, if P = {PI,. . .) PK) = PRIN(A) is the principal component transformation of the albedo time series A, it can be shown that different principal components represent different features in the data. However, it generally cannot be assumed that either the stable components or the temporal components contribute more variance than the other in the original data. The first or the last principal components cannot simply be chosen and used for geological mapping. 3.3. Visual interpretation and separability analysis -

stable component

identi$cation

The principal component analysis offers an engine to separate different components of an image time series. However, identification of which components represent temporal factors and which represent

stable factors needs both visual identification and quantitative analysis. Through visual interpretation, components caused by the presence of clouds, snow, and other significant temporal factors can be identified immediately and excluded from the geological considerations. Since the single image interpretation or analysis may be affected by the presence of clouds and snow in the data, identifying and excluding cloud and snow components provide this method an advantage over the conventional single-shot method. Preliminary selection of the stable components can be visually performed by excluding both temporal and random components. Once an initial selection of the stable components is made, separability analysis can refine the selection of stable components for the classification model. Using the training set data and different variable combinations of the initially selected components, it can be observed which combinations offer the best separability as defined by the Jeffries-Matusita distance. Through this analysis, a subset of the principal components S = {PI,, . . ., Pl#} c (PI, . . ., PK} can be finally selected. Those selected principal com-

D. Yuan et al. /ISPRS Joumul

44

of Photogrammetry & Remote Sensing 53 (1998) 39-53

ponents are more closely related to the lithological background than other components. Maximum likelihood classijication lithological mapping

3.4.

The selected stable components in the study’s experimental classification models were used to perform maximum likelihood classification. Although 18 minimum Euclidean distance, minimum Mahalanobis distance, and maximum likelihood classifiers were tried in the study, extensive visual and statistical accuracy evaluation determined that only the results from the maximum likelihood classifiers were accepted and used for producing the lithological map. The individual techniques or elements used in this model are familiar to most remote sensing scientists. However, this study will demonstrate that those techniques can be organized and used to extract geological information from image time series for lithological mapping.

4. The data set and geological background The complete image data set used in this study consists of 29 MSS images obtained from September 13, 1972, through June 16, 1977. Robinove and Chavez (1978) had studied the albedo temporal variations from a subscene of this MSS series. Through visual evaluation of the quality of the images, and

considering both hardware and software limitations, 14 of the 29 images acquired from January 20 to December 18, 1976, were selected for this geological mapping project. The primary reason for selecting this subscene series was the computer’s limitations in handling a full-scene image time series at the time of the study. Table 1 lists some basic information on the 1976 images used for this study. The selected images covered an area at the border of eastern Nevada and western Utah. The images were georegistered to a USGS base map and subset to cover a selected area centered approximately at West 114”30’and North 30”20’in northern Lincoln County, Nevada. This simultaneously covered area is 300 x 300 pixels (about 23.7 x 23.7 km). Fig. 4 shows the location of the study area and a closer view of the area’s geomorphology. The reference geological map (Fig. 5) was based on the report by Tschanz and Pampeyan (1970). Major geological units are Tertiary volcanic and Quatemary deposit rocks. Table 2 describes these lithological units. The northwest portion of the study area has been largely converted to rangeland and irrigated farmland. Although large portions of the farming and ranching area are in the Quaternary lake deposits (Qr), the straight boundaries of the farmland obviously do not coincide with boundaries between lithologies. This landcover conversion has made the direct interpretation of surface geology from individual MSS scenes difficult. Due to the inaccessibility of the high mountainous area, previous geological

Table 1 MSS image scenes used in this study Sun angle (“)

Visual conditions

76-01-20 76-01-29 76-04- 10 76-05-16 76-06-21 76-07-09 76-08-05 76-08-23 76-09-01 76-10-07 76-10-16

20 24 48 57 58 57 48 44 47 37 30

road, farmland, culture features very clear, snow on mountains road, farmland, culture features very clear, snow on mountains lower left scene smeared, unclear (rain?), snow on mountains vegetation (farmland) area well shown, snow on mountains more blurring on the cultural features whole scene became darker due to haze

76-10-25 76-11-30 76-12-18

32 23 20

Acquisition

date

clear unclear, foggy clear clear jet stream shadow strong, scattered clouds and contract shadows bright scene, vegetation areas well shown bright scene, vegetation areas well shown

D. Yuan et al. /ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

45

Study Aru

ViewUL-- Lake Valley LL-- Fiarview Range

UR-- Wilson Creek Range LR-- Mount Wilson Albedo conversion of MSS data (300x300). Fig. 4. Geographic location and major geomorphic features of the study area. Table 2 Exposed lithological units and descriptions in the study area, where ak are the assigned initial prior probabilities for the maximum likelihood classifiers Age

Symbol

ok

Descriptions

Pliocene Pliocene Tertiary Tertiary Tertiary Tertiary Middle Ordovician Lower Ordovician Uuuer Cambrian

Ql Q0l

0.05 0.50 0.03 0.05 0.01 0.25 0.01 0.01 0.01

Younger lake beds Older alluvium deposits; subdivided into three subunits Qotr, Qor2, and QOrsin the study Younger volcanic rocks, undivided, including intravolcanic sedimentary rocks and perlite Tuffs and tuffaceous sediments Granite stocks and dikes Volcanic rocks, undifferentiated Eureka quartzite Pogonip group Limestone and dolomite

T VY TVt Ts TK,, 0, 0, 6,

mapping for the Mount Wilson region was largely based on air photo interpretation (Tschanz and Pampeyan, 1970). 5. Main processing procedures The data processing procedures can be categorized into three major groups: pre-processing procedures, main processing procedures, and post-pro-

cessing procedures (Fig. 3). We limit this paper’s discussion to the main processing procedure group that consists of multitemporal data reduction, stable component extraction, and geological classification. Obviously, the individual techniques or elements used in the main procedure group are familiar to most remote sensing scientists. However, this study is concerned with the method of organizing those elements to handle multitemporal images and to extract

46

D. Yuan et al./ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

legend 1

Scala

l_-_LY Fig. 5. Base geological map for training the classifier and comparing the results.

geological classifications. It will also show that “the whole is greater than the sum of all components”. Data preparation, such as scene selection, registration, and normalization, can be critical for the success of the whole image time series analysis. Scene selection criteria of multitemporal MSS images for vegetation applications were discussed in detail in Lunetta et al. (1993). One suggestion was to use multiple scenes of the vegetation growing season for vegetation comparison. For geological applications, however, it is better to select scenes of different seasons with the hope that the seasonal effects can be somehow averaged out during the process. Therefore, we used a complete image time series

from different seasons within one year. Multitemporal MSS scene normalization has been discussed by Scott et al. (1988), Hall et al. (1991), Elvidge et al. (1995), and Yuan and Elvidge (1996). In this study, the simplest haze correction technique was adopted for scene normalization. For image registration, a number of control points were digitized from the reference geological map by Tschanz and Pampeyan (1970). All images in the same series were registered to the same map. The root mean square errors (RMSE) of the registration for the image time series are within the limit of one-half (i) pixel (about 40 m). Cubic convolution was used as the resampling technique for the geo-registration of the image series.

D. Yuan et al. /ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

41

5.1. Reduction in spectral domain - albedo conversion

5.2. Component separation in the temporal domain - principal component transformation

Robinove and Chavez (1978) computed the overall average albedos for the images discussed in this paper. The formula for computing albedo for individual pixels is modified from their formula for computing average albedo. As the band contributions to the total or overall albedo have very similar formulas, the band contribution to the total albedo is then redefined as the ‘band albedo’; i.e., the portion of the total albedo that is contributed by electromagnetic radiation in a given band. The albedo conversion formula for the modified band for an individual pixel is: Ci sin(o) where i (i = 1,2,3,4) is the index for a given band, Bi is the band i brightness, Bitin is the minimum value of band i over the whole image, Ci is the conversion factor for band i and the given satellite, a is the sun angle when the image was acquired, and Ai is the albedo value (or, more strictly, the contribution of band i to overall albedo) converted from Bi (Robinove and Chavez, 1978; Robinove, 1982). Actual values for sun angle a! are given in Table 1. Ci actual values are given in Table 3. The overall albedo is a simple sum of all albedo components of all four bands:

Principal component analysis (PCA) is one of the most frequently used data reduction methods for processing multiband images. By comparing the same principal component of one image and that of another, temporal changes in principal component can be detected (Lodwick, 1979). Byrne et al. (1980) superimposed two Landsat images of the same area and treated them as a single eight-channel data set. Ingerbritsen and Lyon (1985) documented the type of behavior observed by Byrne et al. (1980) in two biologically and geologically disparate environments and explained this behavior in terms of the nature of Landsat MSS imagery and the PCA technique. Yuan (1990) investigated the possibilities of decomposing a multitemporal image series into components due to different factors and found that different components usually represent different categories of surface features. For instance, the landscape features are absorbed by the first principal component. Cloud, shadow, and other temporal features are absorbed by the second and the third components. Topographic features are shown in the fourth component, etc. This paper continues the effort of utilizing components for stable geological factors for lithological mapping initiated by Yuan (1990). As a convention,

Aau = AI +

A = (AI, AZ, . . . A14)T

A2

+

A3 +

A4

The corresponding albedo time series of the 1996 MSS time series is shown with appropriate linear stretch in Fig. 6. Seasonal and temporal effects can be easily observed. For instance, the images obtained in winter (76-01-20, 76-01-29, 76-04-10, 76-lo-25,76-ll-30,76-12-18) had snow coverage in the Mount Wilson area. Clouds affected images 7609-10 and 76-10-25. lmage 76-10-16 had a shadow of the jet stream.

Table 3 Constants

used in albedo conversion

for MSS data

Satellites

Cl

c2

c3

c4

Landsat I Landsat II

906.41 881.53

962.03 1131.79

892.6 1 1076.02

392.33 412.98

Source: Robinove

and Chavez,

1978

is denoted as the 14 albedo images obtained through albedo conversion in the spectral reduction procedure. Let P = ET A be the principal component transform for A. It can be shown that Cov(P, P) = A = diag(hl, h2, . . .h14) with h, 5 A2 5 . . . 5 h14. Therefore, it is known that the first PC has the greatest variance, the second PC has the second largest variance, and so forth. The calculated principal component images are shown in Fig. 7. The calculated eigenvalues and their contributions to the total variances are listed in Table 4. Usually the last few PCs have a very small variance, suggesting they represent very little variation in the data set. Those principal components with small variances are conventionally ignored for data reduction purposes. For geological mapping purposes, however, PCs with larger variances may not nec-

48

D. Yuan et al./ISPRS

Journal of Photogrammetry

& Remote Sensing 53 (1998) 39-53

Indices to the albedo images 76-Ol-20 76-06-2 1 76-09-01 76- 1O-25 76-01-29 76-07-09 76-10-07 76-11-30 76-04-10 76-08-05 76-10-16 76-12-18 76-05- 16 76-08-23 Fig. 6. The albedo time series computed from corresponding MSS time series used for this study. Temporal fluctuations due to cloud, snow, and vegetation change can be easily observed.

essarily have the geological information needed, or the needed information may have been smeared by the presence of strong adverse factors. For instance, PC 2 absorbed almost all cloud information from the whole time series with background lithological information suppressed. Similarly snow and the jet stream dominated while the lithological information was suppressed in PCs 3,4, 7, and 8. Although they have relatively larger variance compared with other higher-order PCs, it is intuitively understandable that in the MSS data, the cultural and landscape features, cloud, and snow generate the most variance.

The lithological background, if it can be observed, generates only a small portion of the variance. 5.3. Stable component identijcation - visual and separability analysis The initial stable component identification can be performed by visual inspection. For instance, it is easy to observe that principal components 2, 3,4, 7, 8,9, and 10 are affected to a certain degree by snow, cloud, jet stream shadow, vegetation, etc. Therefore, they should be excluded from consideration. On

D. Yuan et al. /ISPRS Journal of Photogrammetv

& Remote Sensing 53 (1998) 39-53

Indices to the principal PC 5 PC PC ‘1 PC6 PC PC2 PC 7 PC PC3 PC 8 PC4

49

components 9 PC 12 10 PC 13 11 PC 14

Fig. 7. Principal component image series computed from the albedo image series. The first principal component image is basically the averaged image of the albedo series. The second principal component image displays the cloud and shadow effect with background information suppressed.

the contrary, PCs 1, 5, and 6 clearly portray the geomorphology of the area related to the lithological background and are less affected by the temporal factors. The higher-order PCs 11, 12, 13, and 14 contribute only a small portion of the total variance (1.5%). However, certain background information compared to the base geological map still can be observed from the enhanced high-order PCs. For the purposes of experiment, PCs were retained for separability analysis to test if those highorder PCs really help in producing better lithological

classification. Ideal stable components for the classification of the lithology should come from this PC list (1, 56, 11, 12, 13, and 14). Using the lithological training data, the stable components can be further identified together with the separability analysis of the variables. The divergence among samples from different classes can be measured by the Jeffries-Matusita (JM) distance (Swain and Davis, 1978). The JM distance JMii ranges from 0 to 1414. JMij = 0 implies that classes i and j are fundamentally inseparable using the given

50

D. Yuanet al. /ISPRS Journal of Phorogrammetry & Remote Sensing 53 (1998) 39-53

Table 4 The calculated eigenvalues for the albedo time series Eigenvalue

Variance

Variance%

i

4

Ai /d$d

1 2 3 4 5 6 7 8 9 10 11 12 13 14

942.72 180.12 87.40 63.69 47.02 26.35 20.03 13.30 11.71 9.79 7.31 5.69 4.72 3.44

66.23 12.66 6.14 4.47 3.30 1.85 1.41 0.93 0.82 0.69 0.51 0.40 0.33 0.24

Cumulative variance% hk

&I

Aj

/CL:1

kk

66.23 78.99 85.03 89.51 92.81 94.66 96.07 97.00 97.83 98.51 99.03 99.43 99.76 100.00

and maximum likelihood classifiers, were extensively tested using various variables (principal components). The intensive separability analysis, visual evaluation, and the several maximum likelihood classifiers were chosen to produce the composite geological map for the area. Technical details of classic multivariate linear classifiers can be found in many textbooks, such as the one by Schalkoff (1992). The maximum likelihood classifier is given by: @(Pt

4) =

ln(ak) - 0.51n(lVkl) - 0.511X(p, q) - &II&

variables. JMij = 1414 implies that classes i and j are completely separable using the given variables. Table 5 lists the best minimum separability and the best average separability PC combinations from two to six variables out of seven candidates. It was interesting to find out that both PCs 13 and 14 did better than PC 12 in terms of minimum and average separabilities. Therefore, PC 12 was excluded from the stable component list. The remaining list of stable components consisted of PCs 1, 5, 6, 11, 13, and 14.

where ak is the prior probability of class k, X(p, q) is the variable measurement vector at point (p, q), Mk is the mean vector of X within class k, Vk is the covariance matrix of X within class k, IV’1 is the determinant of matrix Vk, and IlX(p, q) - MkllM is the Mahalanobis distance of pixel (p, q) to the center of class k. The Bayesian decision rule for the classification is that if 3b, 1 5 ko < c, where c is the number of classes, such that: D~(P, q) = fini
q)

then pixel (p, q) is classified into class b. The initial estimations of ak values were based on visual interpretations of the base geological map. Table 2 gives these initial estimates. 6. Results and conclusions

5.4. Lithological mapping - maximum likelihood classijcation

In the full-scope study, multivariate linear classifiers, such as minimum Euclidean distance classifiers, minimum Mahalanobis distance classifiers,

The final geological map was a composite from the output of three selected classifiers with combinations 1-5-6-11, 1-5-6-11-13, and 1-5-6-1113-14. The accuracy percentages for those classifiers were all above 85% for the selected training

Table 5 PC combinations with tbe best minimum and average Jeffries-Matusita distance separabilities Number of variables in the combination 2 3 4 5 6

PC combinations with best minimum separability

PC Combinations with best average separability

Best combination

Minimum JM distance

Best combination

Average JM distance

l-6 1-5-11 l-5-6-11 1-5-6-11-13 1-5-6-11-13-14

397 625 736 811 877

l-5 l-5-6 1-5-6-14 1-5-6-11-14 1-5-6-11-12-14

1154 1211 1236 1260 1278

D. Yin

et al. /ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

Fig. 8. The composite

classification

result: a surface geological

data sets. Most classification discrepancies for these classification models occur for class QOl(Quaternary older alluvium deposits). Uncertainties and speckle effects were manually edited to produce the composite. This final map (Fig. 8) can be compared with the base map for classifier training (Fig. 5).

map for the study area.

Since class Q,,t (Quaternary older alluvium deposits) resulted in large discrepancies for these classifiers, three subunits were introduced to resolve the conflicts: Quaternary older lake and alluvium deposits (QOu), Quatemary older alluvium apron deposits (Q,,&, and Quatemary older alluvium fan deposits

52

D. Yuan et al./ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53

(Q&. These conflict areas show significant difference in reflectance in the images, suggesting that these areas differ in terms of material composition, grain sizes, mineral composition, and other physical properties. Therefore it is appropriate to subdivide QOtinto new subunits Q,,u, Qonr and Q,,ts. However, field work has not yet been performed to confirm this subdivision. Multitemporal images can be decomposed into different stable and temporal components for different purposes. The proposed and applied procedure basically consists of four multivariate statistical procedures: the spectral reduction using albedo conversion, component separation using principal component analysis, stable component identification using separability analysis, and lithological classification using the maximum likelihood classifier. When used together, these procedures proved to be effective for extracting and utilizing stable components from multitemporal historical image data. Although the main work presented in this paper was performed five years ago, the importance of using inexpensive historical image data to map geology and other relatively stable features has never been fully recognized in the remote sensing community. While new remote sensing technologies are being developed, the usefulness of existing historical data should not be overlooked. Historical data are often more useful than new data in geological applications, particularly where the surface of the Earth has undergone human and natural disturbances. Historical data are indispensable in environmental change detection. Using historical remote sensing data for geological mapping involves selecting and processing single or multiple scenes. Often tens of MSS scenes are available for an area. The conventional method is to select and work on one scene with the least cloud coverage from multiple choices. But, given the rich information about geology in the historical data, availability of multiple scenes, and the low cost of the data (historical MSS data cost about $200 per scene), it is necessary and possible to take advantage of multitemporal historical image data. Acknowledgements This paper is one in a series of reports of the first author’s work on statistical analysis and inter-

pretation of multitemporal historical Landsat MSS data. Part of this paper was briefly presented at the American Society of Photogrammetry and Remote Sensing annual conference (ASPRS’97) (Yuan et al., 1997). Data used in this paper were provided by Dr. Robinove at the United States Geological Survey. Preliminary work was done under the supervision of Dr. J. Robinson at Syracuse University and Dr. M. Duggin at the State University of New York, College of Forestry and Environmental Sciences (19871992). Manuscript preparation was partially completed when the first author was with the Desert Research Institute (1992-1996). Final completion of the paper was partially supported by the NASA, Commercial Remote Sensing Program Office, under contract number NAS 13-650 at the John C. Stennis Space Center, Mississippi (1996-1997). The authors also wish to thank Marcia Wise and Belle McCann for their valuable technical editing of the manuscript and the NASA CRSP officers for encouraging and permitting us to put this work into a publishable form. Finally, the authors wish to thank the anonymous reviewers of the first draft of this paper for their constructive suggestions and comments that eventually made this paper more suitable for publication. References Byrne, GE, Crapper, EP., Mayo, K.K., 1980. Monitoring landcover change by principal component analysis of multitempoml Landsat data. Remote Sensing Environ. 10, 175-184. Duggin, M.J., 1985. Factors limiting the discrimination and quantification of terrestrial features using remotely sensed radiance. Int. .I. Remote Sensing 6 (l), 3-27. Duggin, M.J., Sakhavat, H., Lindsay, J., 1985. The systematic and random variation of recorded radiance in a Landsat thematic mapper image. Int. J. Remote Sensing 11 (lo), 16691694. Elvidge, C.D., Yuan, D., Weerackoon, R.D., Lunetta, R.S., 1995. Relative radiometric normalization of Landsat multispectral scanner (MSS) data using an automatic scattergram-controlled regression. Photogramm. Eng. Remote Sensing 61 (lo), 12551260. Hall, F.G., Strebel, D.E., Nikeson, J.E., Goetz, S.J., 1991. Radiametric rectification: toward a common radiometric response among multidate, multisensor images. Remote Sensing Environ. 35, 11-27. Ingerbritsen, S.E., Lyon, R.J.P., 1985. Principal components analysis of multitemporal image pairs. Int. J. Remote Sensing 6 (5), 687-696. Justice, CO., Townshend, J.R.G., Holben, B.N., Tucker, C.J.,

D. Yuan et al./ISPRS Journal of Photogrammetry & Remote Sensing 53 (1998) 39-53 1985. Analysis of the phenology of global vegetation using meteorological satellite data. Int. J. Remote Sensing 6, 12711318. Lodwick, G.D., 1979. Measuring ecological changes in multitemporal Landsat data using principal components. Proc. Int. Symp. Remote Sensing of the Environment, Ann Arbor, Michigan, pp. 1131-1141. Lunetta, R.S., Lyon, J.G., Sturdevant, J.A., Dwyer, J.L., Elvidge, C.D., Fens&maker, L.K., Yuan, D., Hoffer, S.R., Weerackoon, R., 1993. North American Landscape Characterization (NALC) Research Plan. EPA/600&93/135, U.S. Environmental Protection Agency. Malingreau, J.P., Belward, A.S., 1992. Scale considerations in vegetation monitoring using AVHRR data. Int. J. Remote Sensing 13 (12), 2289-2307. Pickup, G., Chewings, V.H., Nelson, D.J., 1993. Estimating changes in vegetation cover over time in arid rangelands using Landsat MSS data. Remote Sensing Environ. 43, 243-263. Quarmby, N.A., Milnes, M., Hindle, T.L., Silleos, N., 1993. The use of multi-temporal NDVI measurements from AVHRR data for crop yield estimation and prediction. Int. J. Remote Sensing 14 (2), 199-210. Robinove, C.J., 1982. Computation of physical values for Landsat digital data. Photogramm. Eng. Remote Sensing 48, 781784. Robinove, C.J., Chavez, P.S., Jr., 1978. Landsat albedo monitoring method for an arid region. Presentation AAAS Int. Symp. on Arid Region Plant Resources, Lubbock, Texas. Schalkoff, R., 1992. Pattern Recognition: Statistical, Structural and Neural Approaches. Wiley, New York, p. 364. Scott, J.R., Salvaggio, C., Volchok, W., 1988. Radiometric scene normalization using pseudoinvariant features. Remote Sensing

53

Environ. 26, l-6. Swain, PH., Davis, M., 1978. Remote Sensing: The Quantitative Approach. McGraw-Hill, New York. Tappan, G.G., Tyler, D.J., Wehde, M.E., 1992. Monitoring rangeland dynamics in Senegal with advanced very high resolution radiometer data. Geocarto Int. 1, 87-98. Townshend, J.R.G., Justice, CO., 1986. Analysis of the dynamics of African vegetation using the normalized difference index. Int. J. Remote Sensing 7, 1435-1446. Tschanz, C.M., Pampeyan, E.H., 1970. Geology and mineral deposits of Lincoln County, Nevada. Nevada Bureau Mines Bull. 73, 187. Yuan, D., 1990. Higher order principal components: more information about lithological background. Assn. Petrol. Geochem. Explor. Bull. 6 (l), 49-65. Yuan, D., 1991. Linear temporal models and tests for lithological band brightness in MSS data. Proc. 1991 Int. Geoscience and Remote Sensing Symp., June 3-6, 1991, Helsinki, Finland. Yuan, D., Robinson, J.E., 1990. Time variant effect of vegetation on lithological reflectance. Proc. 1990 Int. Geosci. Remote Sensing Symp, May 20-24, 1990, College Park, Maryland, USA, III: 2113-2116. Yuan, D., Elvidge, C.D., 1996. Comparison of relative radiometric normalization techniques. ISPRS J. Photogranun. 51, 117-126. Yuan, D., Robinson, J.E., Duggin, M.J., 1996. The effect of seasonal factors on geological interpretation of MSS data. Wrster, A., Merriam, D.F. @Is.), Geological Modeling and Mapping, Plenum Press, New York, pp. 153-173. Yuan, D., Holland, D., Lucas, J., 1997. Mapping geology using stable components of multitemporal historical MSS data. Proc. ASPRS ‘97 Convention, April 7-10, 1997, Seattle, WA.