SPECTRAL UNMIXING OF BLENDED REFLECTANCE FOR DENSER TIME-SERIES MAPPING OF WETLANDS

The orbiting cycle and frequent cloud contamination have limited the applications of the moderate-resolution remotely sensed data for detecting rapid land cover changes that are critical to the monitoring of wetlands. It is necessary to use multiple remotely sensed data sources that have different spatial resolution and temporal frequency, because both spatial and temporal details are important in understanding the mechanisms in wetland cover changes. This study examined the applicability of linear spectral mixture analysis to the blended reflectance that was generated by incorporating the enhanced spatial and temporal adaptive reflectance fusion model (ESTARFM). Nine TM and MODIS images of the Poyang Lake area, China acquired in 2004 and 2005 were used to blend the reflectance. In order to account for the spectral variations in materials, we incorporated the multiple endmember spectral mixture analysis (MESMA) in unmixing the blended reflectance. The average absolute differences between the land cover fractions derived from the blended image and those from the observed image were calculated as well as correlation coefficients. Our results demonstrated that MESMA could unmix the blended reflectance generated by ESTARFM. However, due to the existence of the blended pixels with large difference in reflectance from the observed reflectance, the land cover fractions derived from the blended reflectance did not match with those derived from the observed reflectance as well as expected. It is also suggested that the comprehensiveness of the endmember spectral libraries was another factor influencing the agreement.


INTRODUCTION
Taking advantage of regular orbiting intervals and extensive coverage, satellite remote sensing has been utilized as a practical and economical means to monitor and inventory different types of wetlands (Ozesmi and Bauer, 2002). Although a wide variety of time-series remotely sensed data observed with differing sensor designs have been used for the mapping of wetland cover changes, previous studies have shown that wetland mapping using optical remotely sensed data is not as easy as the mapping of other ecosystems (Silva et al., 2008). This is because the spectra of wetland vegetation species show a high level of variability due to the species' structural, biochemical, and biophysical diversity, as well as the spectral confusion among individual wetland components described above (Adam et al., 2010). In addition, due to the tradeoff between spatial resolution and temporal frequency, wetland cover changes has not been monitored with spatial and temporal details simultaneously using the imagery observed by single remotely sensor. For a better understanding of the spatiotemporal dynamics in all land cover components of wetland ecosystems, it is necessary to overcome these difficulties.
In the goal of improving the accuracy of wetland mapping using remotely sensed data, spectral mixture analysis (SMA) have received more attention, due to their relative simplicity of use in * Corresponding author deriving physically interpretable information at subpixel level (Roberts et al., 1993). SMA models mixed spectra in pixels of a remotely sensed image as a combination of endmembers (EMs) − pure spectra representing distinct land cover types . In linear SMA, a spectrum within the instant field of view of a sensor is determined by the sum of each EM spectrum multiplied by its aerial coverage fraction and the residual error. Although many studies have incorporated SMA in the mapping of wetland vegetation and floodplain mapping, only a few studies have also been conducted on SMA using multitemporal remotely sensed data for the mapping of wetland land cover changes. (He et al., 2010;Melendez-Pastor et al., 2010).
In order to increase temporal frequency of moderate-resolution remotely sensed data, several blending techniques have been developed and applied in some studies. Among them, The spatial and temporal adaptive reflectance fusion model (STARFM) (Gao et al., 2006) has been widely used. Recently, Zhu et al. (2010) modified the original STARFM to overcome the poor accuracy of STARFM in heterogeneous landscapes. A few application studies of ESTARFM has proved that the blended reflectance data is comparative to observed reflectance data in chlorophyll index derivation and supervised classification (Singh, 2011;Watts et al., 2011). However, no studies have investigated on the applicability of SMA to the blended reflectance data. In addition, Previous study has not applied these blending techniques in wetland environment.
This study investigates the applicability of SMA to the blended data generated with the enhanced STARFM (ESTARFM) by Zhu et al. (2010) (Zhu et al., 2010) utilizes two pairs of moderateand coarse-resolution data on prior and posterior dates and one coarse-resolution data on the target date. It predicts the surface reflectance of the synthesized moderate-resolution data on the target dates using the linear combination of the spectra for predefined EM land cover classes in the same manner of linear SMA. There are four steps in ESTARFM implementation: (1) Two moderate-resolution scenes are used individually to search for pixels similar to the central pixel in a moving window.
(2) The weights of all similar pixels are determined by the correlation coefficient between moderate-and coarseresolution data (used as a measure of spectral similarity) and geographic distance between the target and similar pixels.
(3) The conversion coefficients are calculated from the surface reflectance of moderate-and coarse-resolution data through linear regression. (4) The surface reflectance of moderate-resolution data on the target date are calculated using the surface reflectance of coarse-resolution data, weights, and conversion coefficients.
Details in the procedure of ESTARFM is described in Zhu et al. (2010).

Multiple endmember spectral mixture analysis (MESMA)
Multiple endmember spectral mixture analysis (MESMA), an extension of SMA, allows EMs to vary on a pixel-by-pixel basis (Roberts et al., 1998). Consequently, MESMA can reduce overall residual error and represent spectral variability in land cover more accurately than conventional linear SMA (Dennison and Roberts, 2003). MESMA is generally implemented by the following procedure: (1) An EM library is constructed from candidate EM spectra.
(2) Optimal EMs are chosen with a EM selection method.
(3) A series of SMA models using user-defined combinations of optimal EMs are applied to every pixel in the image. (4) The model with the minimum root mean square error (RMSE) is selected as the best one from the models that produce physically realistic fractions and meet model conditions. (5) Fractions produced by the optimal models are utilized to map the abundance of EM land cover components. (6) Shade fractions are removed through normalization or addition treatments. (7) EM land cover fractions are validated using higher spatial resolution images or field data. Roberts et al. (2007) describes more details in MESMA implementation.

STUDY AREA
The Poyang Lake area, the largest freshwater lake in China, was selected as the study area in this study ( Figure 1). Poyang Lake (116° 13' E, 29° 9' N) located in the northern part of Jiangxi Province experiences the fluctuation of water level throughout year (Guo et al., 2005). Wetland vegetation in this area is an important food resource for wintering migratory birds, particularly for cranes. It also forms a favorable habitat for Oncomelania snails, the intermediate host for schistosomiasis (Zhou et al., 2005). The dramatic environmental changes in the past decades have consequently made it more difficult to map the change in the distribution of migratory birds and emergence of schistosomiasis. Efficient schemes for its control from the central and provincial government may be difficult to formulate because the effects of the environmental changes in this region (Chen and Lin, 2004). The emergence of highly pathogenic avian influenza, of which migratory birds are believed to be the carrier to poultry birds, is deeply related to the changes in land use and land cover in this region (Feare, 2007).

DATA COLLECTION AND PREPROCESSING
Nine time-series pairs of the Landsat-5 TM Level-1 products and Terra MODIS Daily Reflectance products (MOD09GA) acquired in 2004 and 2005 covering the study area were selected in this study (Table 1). All six spectral bands of TM images except for the thermal band (band 6) and corresponding MODIS bands (bands 1-4, 6, and 7) were utilized. (GCPs) and the nearest neighbor resampling produced an average RMSE of 0.25 pixel. Each of image-to-image registration between the georegistered base image and the other eight TM images produced an average RMSE of less than 0.25 pixel using a first-order polynomial fit with more than 30 GCPs. Brightness values of these geometrically corrected TM images were converted to the ground reflectance through the atmospheric correction using ACORN. The c-correction method (Teillet et al., 1982), a semiempirical function primarily based on cosines of the incident and reflected illumination angles, was applied to minimize the topographic effects. In contrast, only a map reprojection was performed to the time-series MODIS data. This is because MOD09GA products were geometrically accurate when compared with the TM data, and they have already been atmospherically corrected.
Field data of the land cover conditions in the study area were collected in December, 2007 andApril andMay, 2008. In order to comprehensively include all land cover types throughout the study area, the locations of the field record collection were determined with the reference of remotely sensed images in Google Earth. We recorded the GPS coordinates at more than 350 centers of 30 m by 30 m squares covered with single land cover types. Land cover types and existence of land cover changes were also recorded at each location. The field records were classified into four EM land cover classes: (1) green vegetation (GV); (2) non-photosynthetic vegetation (NPV), soils, and impervious surfaces (N/S/I); (3) bright water (W1); and (4) dark water (W2). For each TM image, 240 records consisting of 60 records per EM land cover class were selected from all field records. We chose the field records for image considering the land cover changes by referring the information on land cover changes obtained through the interviews with the local residents and researchers. Therefore, all sets of field records were different for every TM images. These field records for each TM image were utilized as the reference in collecting the candidate image EM spectra and in the class accuracy assessments of the land cover fractions (LCFs) derived through the MESMA of blended reflectance.

Examination of input data combinations for ESTARFM
The agreement between observed reflectance and the reflectance fused with blending techniques is dependent on the combination of input data as Watts et al. (2011) demonstrated. In order to investigate the input data combinations that could achieve highest agreement between the observed and blended reflectance, this research tested all of the possible combinations of the TM and corresponding MODIS images on the TM observation dates. Testing all possible input combinations revealed what data combination was appropriate to blend the relectance in which particular season, particularly when the water level changed rapidly.
The agreement was assessed based on the mean values of the average absolute difference (AAD) between observed and blended reflectance for each band. Since no additional TM data acquired during the studied period was available and the ESTARFM required a prior and a posterior moderate-resolution remotely sensed data, the blended data were generated for the seven TM observation dates: October 28, November 29, and December 15 in 2004, and March 5, August 12, September 13, and September 29 in 2005. The input data combinations followed the chronological order: the prior data were acquired before the target date and the posterior data were acquired after the target date.

Investigation in the applicability of MESMA to blended data
The applicability of MESMA to the blended reflectance was assessed from three perspectives: (1) percentage of modeled pixels; (2) agreement of image dominant land cover classes (EM classes with highest LCFs in pixels) and reference dominant classes; and (3) AAD between the LCFs derived from observed TM data and those from blended data. The blended reflectance generated with the input data combinations that achieved highest agreement between the observed and blended reflectance on the seven TM observation dates were utilized as the blended input data for MESMA.
To establish the approach for EM library development, we made the three comparisons in unmixing the observed blended reflectance using five different EM libraries shown in Table 2. These six cases were set in order to examine whether the candidate EM spectra collected from blended reflectance data need to be included in the candidate spectral library for the MESMA of blended reflectance. Comparisons were made between Case 1 and 2, Case 3 and 4, and Case 5 and 6 utilizing all observed and blended input data on seven target dates.
We utilized the VIPER Tools, a plug-in software under ITT ENVI, in the application of MESMA (Roberts et al., 2007). 40 spectra for each EM land cover class of every images were randomly selected from the image pixels at the locations of 60 field record. The candidate EM libraries were constructed for each input by combining all EM spectra for each case in Table  2. 20 optimal EM spectra (five per land cover class) were chosen from each EM library. using the EM average RMSE (EAR) produced by a spectrum that was used to model all other EM spectra in the same class (Dennison and Roberts, 2003). one spectrum of photogrammetric shade was added to account for the spectral variation of reflectance data. When a series of SMA models were applied to every pixel of each reflectance data, two-, three, and four-EM models were applied in this study. The combinations of EM land cover classes used in this study are shown in Table 3. Maximum RMSE (2.5% reflectance), maximum fraction (105% reflectance), minimum fraction (-5% reflectance), and maximum shade fraction (50% reflectance) restrictions were applied to the SMA models (Michishita et al., 2012). In the selection of optimal models for every pixels, we empirically determined the RMSE threshold between two-and three-EM models as 0.2%, that between three-and four-EM models as 0.4%, and that between two-and

Input data combinations for ESTARFM
The AADs between the observed and blended reflectance for the input combinations that achieved highest mean values of AADs in all spectral bands on each target date are shown in Table 4. The three alphabet codes stand for the input combinations. Each alphabet in the codes corresponds to observation dates as summarized in Table 1. The first and third lower alphabets refer to the prior and posterior dates, and second upper alphabets refer to the target dates. Therefore, one TM and one MODIS images were utilized on the dates denoted by lower case letters, while one MODIS image was utilized on the dates denoted by upper case letters.
The tests of all input combinations demonstrated that ESTARFM did not always obtain strong agreement when the prior and posterior dates are closer to the target dates. This fact suggested that the correlation between the observed and blended reflectance may have depended on other environmental factors in the study area, such as phenological stage of the wetland vegetation and water level of the lake. The difference between observed and blended reflectance was generally larger in the spectral bands with longer wavelengths, such as nearinfrared (NIR), shorter and longer shortwave infrared (SWIR) bands. Larger differences were found particularly in the visible blue and NIR bands when the input reflectance with cloud contamination was utilized to blend the reflectance.
6.2 Applicability of MESMA to blended data 6.2.1 Percentages of modeled pixels: The percentages of pixels modeled by MESMA for the observed and blended reflectance data is shown in Figure 2. All three cases using the observed reflectance data (Cases 1, 3, and 5) achieved high percentages of modeled pixels (higher than 89%) except for Case 3 on September 13, 2005. However, many pixels in the blended data could not be unmixed by other three cases (Cases 2, 4, 6), particularly on August 12 and September 13 in 2005.  Figure 2. Percentages of modeled pixels in the MESMA of observed and blended reflectance data This is mainly because the input blended data for those dates had larger differences in reflectance from the observed data in the blue, green, read and NIR bands than the blended reflectance for other dates. It led to the RMSEs larger than the modeling threshold (2.5% reflectance) in all EM models. For each comparison, the cases for the observed data reached higher percentage of modeled pixels than those for the blended data, except for in the Comparison 1 on March 5, 2005 and Comparison 2 on September 29, 2005. Although no significant difference in the percentage of modeled pixels was found among the three cases for blended data, differences in the percentage of modeled pixels (larger than 5%) were confirmed among the three cases for observed data on August 12 and September 13 in 2005, suggesting that the differences were brought by the inclusion of EM spectra on different dates and data sources (observed, blended or both reflectance). Figure 3. Kappa coefficients in EM class accuracy assessment 6.2.2 Class accuracy assessment: Figure 3 illustrates the kappa coefficients in the accuracy assessment of the dominant EM classes in the LCF images derived from the observed and blended reflectance data with different EM libraries. Kappa coefficients were calculated in the same manner as in the accuracy assessment of hard-labeled classification maps using confusion matrices. Low kappa coefficients were chiefly brought by the unmodeled pixels, unmixing of the blended reflectance with unrealistic EM combinations, confusion between two water classes, and overestimation and underestimation of the LCFs from the blended reflectance due to the overestimation and under-estimation of reflectance values. Although Case 5 (observed EM library for observed data) achieved kappa coefficients higher than 0.86 on all seven dates, kappa coefficients for other five case varied on different dates. Kappa coefficients for the LCFs derived from observed data was generally higher than those derived from blended data in Comparisons 1 and 3. Low kappa coefficients for Case 1 on September 13, 2005 and those for Cases 2, 3, 4, and 6 on March 5, August 12, September 13, and September 29 in 2005 implied the difficulty in collecting the candidate image EM spectra comprehensively in these seasons. Case 3 had lower kappa coefficients than Case 1 on all dates except for October 28, 2004, indicating that the inclusion of too many candidate observed EM spectra on different dates that had true physical values led to the reduction in the representativeness of the optimal EM spectra. In contrast, kappa coefficients for Cases 4 and 6 were as similar as or higher than those for Case 2 on most dates, implying that the observed and blended data had large difference in their spectral characteristics.

Average absolute difference:
The AADs between the LCFs derived from the observed and blended data were summarized in Table 5. Pixels which were unmodeled in either case of the comparisons were excluded from the calculation. W2 had the highest AADs of LCFs on all seven dates in the three comparisons, obtaining higher than 22% of fractions. Table 5. Average absolute differences between the LCFs derived from the observed and blended reflectance data Mean values of the AADs for each EM class over time revealed the large difference (larger than 11%) between LCFs derived from observed reflectance and those derived from blended reflectance. The temporal variations of the AADs for the two ground EM land cover types (GV and N/S/I) were smaller than those for the two water land cover types (W1 and W2) in all comparisons. In Comparisons 2 and 3, the mean values of the AADs over the land cover classes on the four dates in 2005 were smaller than those on the three days in 2004, suggesting that the optimal EM spectra chosen from blended candidate spectra were not so representative as those chosen from observed candidate spectra. When mean values of the AADs of LCFs for all EM land cover classes were compared, Comparison 1 using the EM spectra collected from the observed data on the prior and posterior dates achieved the highest values of the three comparisons. In summary, the LCFs derived from the blended data achieved strongest agreement with those derived from the observed reflectance data when the optimal EM spectra were chosen from the candidate spectra collected from the observed reflectance data on prior and posterior dates of the target dates.

CONCLUSIONS
This study investigated the applicability of MESMA to the blended reflectance data generated with ESTARFM utilizing Although ESTARFM could blend the reflectance accurately in visible and shorter SWIR bands, it could not always generate accurate blended reflectance in NIR and longer SWIR bands. The large differences between the observed and blended reflectance in those bands brought a large number of unmodeled pixels, unmixing of the blended reflectance with unrealistic EM combinations, confusion between two water classes, and overestimation and underestimation of the LCFs. And they consequently led to the poor percentages of modeled pixels, low kappa coefficients in EM class accuracy assessment, and large AADs between the LCFs derived from the observed and blended reflectance. Although there is a need for the refinement in building the EM library, this study achieved strongest agreement between the LCF derived from the observed and blended reflectance data when using the optimal EM spectra chosen from the candidate spectra collected from the observed data on prior and posterior dates of the target dates.