Next Article in Journal
Hub Height Ocean Winds over the North Sea Observed by the NORSEWInD Lidar Array: Measuring Techniques, Quality Control and Data Management
Previous Article in Journal
Recent Declines in Warming and Vegetation Greening Trends over Pan-Arctic Tundra
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Remote Sensing
Volume 5
Issue 9
10.3390/rs5094255
remotesensing-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Mapping and Evaluation of NDVI Trends from Synthetic Time Series Obtained by Blending Landsat and MODIS Data around a Coalfield on the Loess Plateau

by
Feng Tian
1,2,
Yunjia Wang
2,*,
Rasmus Fensholt
3,
Kun Wang
4,
Li Zhang
4 and
Yi Huang
2
1
Jiangsu Key Laboratory of Resources and Environmental Information Engineering, China University of Mining and Technology, No.1 Daxue Road, Xuzhou 221116, China
2
School of Environment Science and Spatial Informatics, China University of Mining and Technology, No.1 Daxue Road, Xuzhou 221116, China
3
Department of Geosciences and Natural Resource Management, University of Copenhagen, Øster Voldgade 10, DK-1350 Copenhagen, Denmark
4
Key Laboratory of Digital Earth Science, Institute of Remote Sensing and Digital Earth, Chinese Academy of Sciences, No.9 Dengzhuang South Road, Beijing 100094, China
*
Author to whom correspondence should be addressed.
Remote Sens. 2013, 5(9), 4255-4279; https://doi.org/10.3390/rs5094255
Submission received: 24 June 2013 / Revised: 24 August 2013 / Accepted: 26 August 2013 / Published: 2 September 2013
Browse Figures
Versions Notes

Abstract

: The increasingly intensive and extensive coal mining activities on the Loess Plateau pose a threat to the fragile local ecosystems. Quantifying the effects of coal mining activities on environmental conditions is of great interest for restoring and managing the local ecosystems and resources. This paper generates dense NDVI (Normalized Difference Vegetation Index) time series between 2000 and 2011 at a spatial resolution of 30 m by blending Landsat and MODIS (Moderate Resolution Imaging Spectroradiometer) data using the Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM) and further evaluates its capability for mapping vegetation trends around a typical coalfield on the Loss Plateau. Synthetic NDVI images were generated using (1) STARFM-generated NIR (near infrared) and red band reflectance data (scheme 1) and (2) Landsat and MODIS NDVI images directly as inputs for STARFM (scheme 2). By comparing the synthetic NDVI images with the corresponding Landsat NDVI, we found that scheme 2 consistently generated better results (0.70 < R2 < 0.76) than scheme 1 (0.56 < R2 < 0.70) in this study area. Trend analysis was then performed with the synthetic dense NDVI time series and the annual maximum NDVI (NDVImax) time series. The accuracy of these trends was evaluated by comparing to those from the corresponding MODIS time series, and it was concluded that both the trends from synthetic/MODIS NDVI dense time series and synthetic/MODIS NDVImax time series (2000–2011) were highly consistent. Compared to trends from MODIS time series, trends from synthetic time series are better able to capture fine scale vegetation changes. STARFM-generated synthetic NDVI time series could be used to quantify the effects of mining activities on vegetation, but the test areas should be selected with caution, as the trends derived from synthetic and MODIS time series may be significantly different in some areas.

1. Introduction

With an area of 648,700 km2, the Chinese Loess Plateau is home to about 108 million people [1]. It is one of the largest regions with severe soil erosion in the world [2]. The largest quantity of coal resources in China are produced in the Loess Plateau (Figure 1), and the total output is increasing. Coal mining activities (coal exploiting) accelerate soil erosion, due to the alteration of underground geological structures and land surface subsidence [3,4]. It is important to assess the effects of coal mining activities on the local environment. One solution to this issue is to compare the changes in vegetation productivity, which is recognized as one of the key indicators of environmental conditions, between mined and non-mined areas [5]. Therefore, long-term earth observation data with high spatial and high temporal resolution are of great interest for the assessment of environmental conditions in coalfields on the Loess Plateau.

Landsat TM/ETM+ (Thematic Mapper/Enhanced Thematic Mapper Plus) data with a 30-m spatial resolution seems to be the best data to capture vegetation complexity at a local scale compared to 4-m IKONOS (Space Imaging) data and 250-m MODIS (Moderate Resolution Imaging Spectroradiometer) data [6] and have been widely used to map vegetation changes over local and regional scales [7–12]. Landsat scenes from the same month or season of each year are used to assemble time series data for vegetation trend analysis [8,9]. However, this kind of time series often suffers from missing data in a certain year and false trend results that are caused by variation in observation time, due to the 16-day Landsat revisit cycle and frequent cloud contamination. In addition, Landsat data alone cannot be used for phenology analysis, because of the sparse and unbalanced distribution of the acquisition date [8], limiting its application in monitoring of long-term vegetation phenology change.

In contrast, satellite sensors, such as the Advanced Very High Resolution Radiometer (AVHRR), Système Pour l’Observation de la Terre (SPOT)-VEGETATION, MODIS and ENVISAT-Medium Resolution Imaging Spectrometer (MERIS), are able to acquire observation data more frequently. Cloud contamination can be reduced or even removed by compositing images within a certain number of continuous days, thus allowing robust trend and frequency analysis to obtain phenological information [8]. With these more frequent time series data, many vegetation trend analyses have been carried out at regional to global scales [13–20]. However, the spatial resolution (250–8,000 m) of these remote sensing data is too coarse for vegetation assessment at local or smaller scales [8,21].

With a spatial resolution of 250 m or 500 m, the MODIS images can be used for monitoring the overall changes of vegetation in a whole mining area. However, the width of a single coal face is usually 100–300 m and the length, 500–2,000 m. At such a scale, Landsat images with a spatial resolution of 30 m are suitable to capture vegetation changes, while the MODIS pixel size is too coarse. Generally, the coal cutting speed is 1–2 m per day (30–60 m per month); consequently, vegetation change monitoring at a temporal resolution better than one month is necessary. Therefore, it is of great interest to combine the spatial resolution of Landsat images and the temporal resolution of MODIS images. In addition, improvements in resolution are also highlighted in mapping land surface temperature [22–24].

Several data fusion methods have been proposed to gap-fill the trade-off of remote sensing sensors between spatial and temporal resolution [21,25–28], among which the Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM), proposed by Gao et al. [21], has been recognized as an efficient one [29]. The feasibility of using STARFM to blend Landsat and MODIS surface reflectance data to generate synthetic reflectance data with high spatial and high temporal resolution has been shown by several studies in various regions with different vegetation types [30–35]. Several studies modified the STARFM algorithm to better capture spatial and temporal land cover changes [28,36–38]. Singh [31] generated a long-term Normalized Difference Vegetation Index (NDVI) time series of cropland in India and found that the synthetic time series adequately captured the phenological patterns over a growing season in eight separate years. Bhandari et al. [33] and Schmidt et al. [35] evaluated the ability of STARFM to generate long-term synthetic time series data to extract phenological parameters using TIMESAT (a software for analyzing time-series of satellite sensor data) [39] and obtained reasonably accurate results for different vegetation communities in Australia. However, few applications of long-term vegetation trend analysis have been carried out using the STARFM-generated synthetic NDVI time series data.

This study is expected to generate synthetic NDVI time series at a spatial resolution of 30 m using STARFM and to evaluate its capability for monitoring vegetation changes around a typical coalfield on the Loess Plateau. Landsat and MODIS NDVI images were used in a two-fold way: (1) directly as inputs for STARFM to generate synthetic 30-m resolution NDVI time series and (2) for evaluation of the accuracy of the synthetic NDVI time series produced. Synthetic 30-m spatial resolution NDVI images were evaluated against four independent Landsat scenes covering different growth stages in 2000. Time series data from 2000 to 2011 were generated based on a limited number of Landsat scenes, and trend analyses were conducted for both the synthetic dense NDVI time series and annual maximum NDVI (NDVImax) time series. The accuracy of these high spatial resolution trends was evaluated by comparing with the corresponding trends derived from MODIS NDVI time series. The methods presented provide a framework for the analysis of trends from synthesized data using the STARFM algorithm.

2. Data and Methods

2.1. Study Area

The Datong Coalfield is one of the largest coalfields in China, which is located in Shanxi Province in the eastern part of the Loess Plateau (Figure 1). The geomorphic landscape in the western part of this region is characterized by ravines and gullies, which are typical on the Loess Plateau. In contrast, the landscape in the eastern part, where Datong city is located, is much flatter. The main vegetation types are irrigated and non-irrigated croplands, shrubs, broadleaf forest and grassland. The land cover map (Figure 2) was created using the Support Vector Machine (SVM) method by Zhao et al. [40], having an overall accuracy of over 85%. According to this classification map, the accuracy of STARFM-generated synthetic data will be assessed among different vegetation types. For this purpose, the 85% classification accuracy is adequate. The growing season of the natural vegetation types lasts from April to November, while farming periods are much shorter, from June to August. The annual precipitation in this region is about 250–450 mm.

Coal mining activities in the Datong Coalfield has a long history, dating from ancient China. Mining occurs in the Jurassic and Carboniferous-Permian strata, with areas of 772 km2 and 1,739 km2, respectively. Until recently, all the mining activities were located in the Jurassic strata, as shown by the boundary in Figure 1. The majority of coal resources in these strata were found about 300 m below the surface, and about 690 km2 of land has experienced surface subsidence, due to the mining activity. As a result of this subsidence and subsequent soil erosion, trees cannot live in these areas.

2.2. MODIS Data

The MODIS vegetation index product, MOD13Q1, with a 16-day composite and a spatial resolution of 250 m, was selected to provide temporally frequent coarse resolution data. It contains both the NDVI product and reflectance layers of near infrared (NIR) and red that will be used in the following analysis. All MOD13Q1 scenes of h26v04 and h26v05 from 2000 to 2011 were obtained from the National Aeronautics and Space Administration (NASA) Earth Observing System Data and Information System ( http://earthdata.nasa.gov). To meet the STARFM input requirements, MOD13Q1 data was subset to the extent of the study area, re-projected to the same Universal Transverse Mercator (UTM) projection as Landsat and resampled to a resolution of 30 m using the MODIS re-projection tool (MRT). A nearest neighbor resample approach was used to ensure the least distortion of the original NDVI value.

The MOD13Q1 vegetation product contains a summary quality layer called “pixel reliability”, holding information on the overall pixel quality. In this study, only pixels with flag values equal to 0 (good data) were used, whereas pixels with other flag values were masked out. A detailed explanation of MOD13Q1 pixel reliability can be found at http://www.ctahr.hawaii.edu/grem/mod13ug/sect0005.html. As MODIS data were not collected in the first two months of the year 2000, average values of the same period between 2001 and 2011 were used to form a completed time series that was needed for the trend analysis [41].

A simple five-point temporal filter was applied to the MOD13Q1 NDVI time series to fill gaps and missing values from bad quality pixels in the original time series [41–43]. For a missing NDVI value at position X in the time series, the average of NDVI values at X − 1 and X + 1 is given. If only one of the two values is available, then this value will be used. If both the values are missing, NDVI values at X − 2 and X + 2 will be included in the filtering process, following the same rules as for X − 1 and X + 1.

2.3. Landsat Data

Landsat TM/ETM+ images L1T (Level 1 Terrain Corrected) products of WRS (Worldwide Reference System)-2, Path 125/Row 32 were used. A total number of 32 Landsat scenes were available (cloud cover <10% and SLC (Scan Line Corrector) -on) from 2000 to 2011 and all of them were downloaded using the United States Geological Survey (USGS) GLOVIS (Global Visualization Viewer) portal ( http://glovis.usgs.gov). These Landsat images were found to be spatially consistent at a high level in the study area, thus no additional registration work was done. Atmosphere correction of Landsat images was performed using the LEDAPS (Landsat Ecosystem Disturbance Adaptive Processing System) tool [44], in which column water vapor is taken from NOAA (National Oceanic and Atmospheric Administration) National Centers for Environmental Prediction (NCEP) reanalysis data, ozone concentrations are derived from Total Ozone Mapping Spectrometer (TOMS) data, and aerosol optical thickness (AOT) is directly extracted from the imagery using the Dark Dense Vegetation (DDV) method [45]. As the LEDAPS employs the same radiative transfer model for MODIS atmospheric correction, the difference between Landsat and MODIS reflectance data was minimized. The NIR and red reflectance bands in Landsat were then extracted and subset to the extent of the study area. Under ideal conditions, atmospheric-corrected Landsat observations should be identical to the MODIS surface reflectance products, which is one of the assumptions of the STARFM algorithm. However, systematic deviation between the two datasets will always exist [31,33] that could affect the performance of STARFM. To minimize this impact, Landsat NIR and red reflectance bands were normalized to the corresponding MODIS images by subtracting the mean difference value between them.

2.4. Generation of Synthetic NDVI Time Series Using STARFM

2.4.1. STARFM

Neglecting geolocation errors and differences in atmospheric correction, a coarse-resolution MODIS pixel can be aggregated from the corresponding fine-resolution Landsat pixels. For a homogenous pixel at the MODIS resolution, Landsat surface reflectance can be obtained by adding an error value (caused by differing bandwidth and solar geometry) to the MODIS surface reflectance [21]. Suppose the land cover type and system errors do not change over time; we will be able to predict the Landsat surface reflectance on a specific day (tn) by blending one or more pairs of Landsat and MODIS surface reflectance images on a base day (t0) and the MODIS surface reflectance image acquired on day tn [21]. However, the relationships between MODIS and Landsat reflectance are complicated by the fact that MODIS observations may be mixed pixels, and land cover type and Bidirectional Reflectance Distribution Function (BRDF) may change over time [21]. To solve this problem, STARFM employed a moving window to incorporate additional information from neighboring Landsat and MODIS image pixels using a weighting function defined by a combination of spectral, temporal and spatial differences among the given Landsat and MODIS scenes [21]. Surface reflectance time series can be obtained at the same spatial resolution as Landsat and the same temporal resolution as MODIS using the STARFM algorithm. The synthetic work was performed using the STARFM tool acquired from the Landsat Ecosystem Disturbance Adaptive Processing System ( http://ledapsweb.nascom.nasa.gov/).

2.4.2. Two Schemes of using STRFM

Two schemes were proposed to generate a synthetic time series of 30-m resolution NDVI images using STARFM (Figure 3). Scheme 1: The NIR and red reflectance bands of Landsat and MODIS images on day t0 and MODIS bands (NIR and red) on day tn were used as inputs for STARFM to generate synthetic bands (NIR and red) on day tn. Then, the NDVI image on day tn was calculated using the synthetic NIR and red bands generated above. Scheme 2: the Landsat NDVI image on day t0 was first derived from NIR and red bands, which was then used directly as an input for STARFM, together with MODIS (MOD13Q1) NDVI images on day t0 and tn, to generate the synthetic NDVI image on day tn. Both of these two schemes were used to generate synthetic NDVI images, and the accuracy of the results were evaluated by comparing with Landsat NDVI on day tn. A random sample (25%) of the vegetation pixels in both independent Landsat and synthetic NDVI images were compared on a per-pixel basis by means of correlation. The scheme with higher accuracy was chosen to generate synthetic NDVI time series from 2000 to 2011 and was further used for the following trend analysis.

2.4.3. Selection of Base Landsat Images

Table 1 shows the base pair images used in STARFM to generate the synthetic NDVI time series from 2000 to 2011. Although 32 Landsat scenes were available in this study area, nearly half of them were obtained in the non-growing season, and the temporal distribution of these scenes was very unbalanced. Consequently, one base pair of Landsat and MODIS images were selected for each year. For those years where no Landsat images were available, base pair images from the adjacent year were used.

2.5. Accuracy Assessment of STARFM Synthetic Data

Available Landsat images in 2000 were used for accuracy assessment of the two schemes above, as there were five cloud-free Landsat scenes covering the whole growing season (6 May 2000, 01 July 2000, 11 September 2000, 13 October 2000 and 21 October 2000). MODIS images with a composite period closest to the observation date of Landsat scenes were selected to form image pairs. The Landsat image of 13 October 2000 and the MODIS composite image from 29 September 2000 to 14 October 2000 were used as the base pair images for STARFM to generate synthetic images corresponding to the other four dates (6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000).

To evaluate the accuracy of synthetic NDVI time series, the mean and standard deviation values of synthetic dense NDVI time series images and NDVImax time series images were calculated and compared with the corresponding MODIS time series.

2.6. Trend Analysis and Evaluation

Trend analysis of NDVI time series was carried out in two stages. First, all of the synthetic NDVI images were included in the analysis to ensure that the temporal properties of the synthetic time series data were fully utilized. Theil-Sen (TS) slope trend analysis was used to quantify the temporal trends in the NDVI time series. This is a non-parametric statistical method calculating the median of the slopes between all n ( n − 1 ) 2 pair-wise combinations over time [46,47]. The TS slope can tolerate up to approximately 30% noisy input observations without influencing the trend estimate [48] and is robust against seasonality, non-normality and serial dependence [49,50]. The significance of the NDVI time series trend was examined by the non-parametric seasonal Mann-Kendall (SMK) test, which is a completely rank-based test that is robust against non-normality, missing values, seasonality and, if corrected, can be robust against serial dependence [49,51,52]. The SMK test produces Z scores to measure the significance of a monotonic trend over time, with a Z value >1.96 representing a significant increasing trend and a Z value <−1.96 indicating a significant decreasing trend (p < 0.05). In this study, SMK Z scores were obtained using the method proposed by Hirsch et al. [49], considering serial dependence in the seasonal time series. This method was recommended by de Beurs et al. [51] for remote sensing time series analysis and has been used for vegetation trend analysis with dense time series satellite data [53–55].

In the second stage, the annual maximum NDVI (NDVImax) values for every pixel in each year between 2000 and 2011 were extracted from both synthetic and MODIS time series data and used for trend analysis. NDVImax values represent the highest vegetation productivity in one year, thus avoiding seasonal variation and minimizing the errors introduced by STARFM. In this stage, temporal trends in the NDVImax time series were also examined by the TS slope. The Mann-Kendall (MK) test [56,57], which does not consider seasonal variation, was used to calculate the significance of the trends.

Both the NDVI and NDVImax trends derived from STARFM synthetic time series were evaluated by comparing to the corresponding trends derived from the MODIS time series.

3. Results

3.1. Synthetic NDVI Time Series

Figure 4 shows the scatterplots between Landsat NDVI images of 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000 and the corresponding synthetic NDVI images. Figure 4A shows the results of scheme 1 using band (NIR and red) reflectance images on 13 October 2000 as inputs for STARFM, while Figure 4B shows the results of scheme 2 using NDVI images on 13 October 2000 directly as an input for STARFM. For the four validation images in this study area, the results of scheme 2 showed a higher accuracy, with scatterplots closer to the 1:1 line and R2 > 0.70, in contrast to those of scheme 1, in which 0.56 < R2 < 0.70. The figures for different vegetation types in Table 2 give similar results, and scheme 2 consistently generates higher R2 values and lower mean absolute difference values than scheme 1. Consequently, scheme 2 was selected to generate synthetic NDVI time series from 2000 to 2011.

To clarify whether using base pair images from the adjacent year would introduce significant errors, base pair Landsat and MODIS images on 11 October 2002 were used as inputs for STARFM to generate synthetic NDVI on 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000. For the four validation Landsat NDVI images, scatterplots and R2 values of STARFM results using base pair images two years later (11 October 2002) were comparable to those using base pair images in the same year (13 October 2000) (Figure 4B,C). Similar results are shown in Table 2 for R2 values and mean absolute difference for different vegetation types. Therefore, with base pair images from the adjacent year, STARFM could still maintain the accuracy of synthetic NDVI images.

Synthetic and MODIS NDVI images were compared in terms of the mean and standard deviation of dense time series images and NDVImax time series images during the period of 2000–2011 (Figure 5). The mean values in the two time series showed a high level of consistency, while synthetic NDVI images show slightly higher standard deviation values than MODIS images. It is likely that synthetic images maintained the high spatial resolution characteristics from Landsat images and captured, at the same time, NDVI changes in the MODIS time series. The NDVImax images using the same base Landsat image (Table 1) showed similar standard deviation values (Figure 5B), possibly since the spatial patterns of the synthetic images were determined by the base Landsat image [21].

3.2. Comparison of Trends from Synthetic and MODIS Time Series

The trends from NDVI dense time series and annual NDVImax time series during the period of 2000–2011 were obtained from both synthetic and MODIS (MOD13Q1) time series at a spatial resolution of 30 m and 250 m, respectively. Using MODIS as the reference time series, comparisons between the two results were carried out, to evaluate the accuracy of the NDVI trends derived from the synthetic time series.

3.2.1. Results of Dense NDVI Times Series

Figure 6 shows the spatial distribution of NDVI trends derived from both synthetic and MODIS dense time series for the entire study area during the period of 2000–2011. Figure 6A,B shows the magnitude of NDVI trends over time (NDVI slopes shown in NDVI units per year), while Figure 6C,D shows the statistical significance of the trends, according to the SMK test (Z scores). The results of the synthetic time series shown in Figure 6A,C exhibit a similar overall pattern with the results of the MODIS time series presented in Figure 6B,D. With a spatial resolution of 30 m, however, the results of the synthetic time series provide more detailed information than the results of the MODIS time series, with a spatial resolution of 250 m. More areas with a negative trend were found in the results of the synthetic time series than in the results of the MODIS time series, especially in the mined area, although positive trends dominated for both time series.

Figure 7 shows the scatterplots of the SMK Z scores versus NDVI slopes for the vegetated pixels in the study area, indicating the relationship between the magnitude of NDVI trends and the statistical significance of the trends. The total number of pixels in the NDVI trends derived from the synthetic time series is nearly 70-times larger than in those derived from the MODIS time series, which equals the ratio of per pixel area between MODIS and Landsat images. In both scatterplots, a small proportion of pixels have a positive slope, but a negative Z value (between −1.96 and zero in the second quadrant). The discrepancy between TS and SMK was also reported by Fensholt et al. [54], who used a prewhitening method [58] to correct the serial correlation within the seasonal time series. However, this inconsistency in the results of the two statistical methods could be tolerated, because none of the positive TS slopes in the second quadrant passed the SMK significance test. This explains the reason for the more positive trends in the NDVI slope maps (Figure 6A,B) than in the corresponding significance maps (Figure 6C,D). The majority of the high frequency pixels in Figure 7A show a similar pattern with those in Figure 7B in terms of the relationship between Z scores and TS slopes, while the NDVI slopes in Figure 7A extend to wider ranges. This may be caused by the higher spatial resolution of the synthetic images, which allows them to capture NDVI changes at the scale smaller than the spatial resolution of the MODIS images.

The frequency distributions of NDVI slopes in Figure 6A,B were further studied and are shown in Figure 8. The mean values of the NDVI slopes of both NDVI trends were nearly the same, with 0.0036 for those derived from synthetic time series and 0.0031 for those derived from MODIS time series. However, the standard deviation value was much higher for NDVI slopes derived from synthetic time series (0.0057) than those derived from MODIS time series (0.0020). Due to the huge difference of the standard deviation values and the total pixel numbers between these two distributions, both the means and standard deviation values of them were found to be significantly (p < 0.001) different from each other according to the t-test (two-side) and F-test, respectively. In both histograms, pixels characterized with slope values between 0.002–0.004 account for the largest proportion. However, the percentages of these pixels differed substantially, with 27.6% of the total vegetated pixels in the NDVI slopes concentrated between 0.002–0.004 derived from the synthetic time series, while 63.1% in that derived from the MODIS time series. Possible reasons for these differences are discussed in Section 4.2.

Table 3 reflects similar spatial patterns in the two NDVI trends, showing a cross-tabulation with the vegetated surface area characterized by positive and negative trends, according to the SMK test. The p-value obtained from a chi-square test of the cross-tabulation table was less than 0.001, meaning that there was a statistically significant relationship between the two trends [59]. The contingency coefficient, which is recognized as an objective measure for the degree of dependence of the classifications in the cross-tabulation table, was 0.56. For a 4 × 4 matrix, the maximum contingency coefficient is 0.86. Consequently, the NDVI trends derived from synthetic dense time series showed a high level of consistency with the NDVI trends derived from MODIS time series, in terms of the significance of the NDVI trends obtained by the SMK test.

3.2.2. Results of NDVImax Time Series

Figure 9 shows the spatial distribution of the NDVImax trends derived from both synthetic and MODIS NDVImax time series during the period of 2000–2011. Figure 9A,B shows the magnitude of NDVImax trends over time (NDVImax slopes in NDVI units per year), while Figure 9C,D shows the statistical significance of the trends, according to the MK test (Z scores). The patterns between the two NDVImax trends were visually consistent at a high level.

Figure 10 shows the scatterplots of the MK Z scores versus NDVImax slopes for the vegetated pixels in the study area, and similarities in the patterns of the relationship between Z scores and NDVImax slopes were evident. All of the pixels are located in the first and third quadrants, showing the consistent results of the TS slope estimator and the MK significance test. The frequency distributions of NDVImax slopes in Figure 9A,B are shown in Figure 11. The mean values of NDVImax slopes of both trends were nearly the same, with 0.0009 for those derived from synthetic time series and 0.0007 for those derived from MODIS time series. The standard deviation was higher for NDVImax slopes derived from synthetic time series (0.0090) than for those derived from MODIS time series (0.0065). Like the results of Figure 8, both the means and standard deviation of these two distributions were significantly (p < 0.001) different according to the t-test (two-side) and F-test, respectively. The possible reasons for these differences are discussed in Section 4.2.

The cross-tabulation of the results of the NDVImax trends is shown in Table 4, illustrating the vegetated surface area characterized by positive and negative trends, according to the MK test. A statistically significant relationship was found between the two NDVImax trends, with a p-value less than 0.001 and a contingency coefficient of 0.57 out of 0.86 obtained from the chi-square test of the cross-tabulation table. Therefore, the NDVImax trends derived from synthetic time series were highly consistent with the NDVImax trends derived from the MODIS time series, in terms of the significance of NDVImax trends obtained by the MK test.

4. Discussion

4.1. Accuracy of the Synthetic NDVI Images

The accuracy of STARFM synthetic NDVI images derived from scheme 2, in which NDVI images were used directly as input, was higher than that derived from scheme 1, in which band reflectance images were used as input (Figure 4A,B). This may be because the nonlinear relationships between NDVI and reflectance bands amplified the errors in the synthetic NIR and red reflectance images when transforming into NDVI images, which was also reported in the results of Walker et al. [32]. However, the decrease of synthetic accuracy in NDVI was seldom seen in the results of Bhandari et al. [33] in terms of the R2 values, indicating that the performance of the two schemes may differ between separate study areas. The huge numbers of ravines and gullies in the study area increase the variation of BRDF (Bidirectional Reflectance Distribution Function) and the heterogeneity of the land surface, which may explain the lower accuracy of synthetic NDVI images (Figure 4, Table 2) compared to other studies, where the landscapes were more flat and homogeneous [30,32,33]. Bhandari et al. [60] examined the impact of solar and viewing geometry on the MOD13Q1 product in an Australian forest and found that there existed a higher amplitude and a phase shift in the MOD13Q1 NDVI time series compared to NDVI derived from the BRDF adjusted MODIS product (MCD43A4). Hence, using the MOD13Q1 product as input for STARFM may introduce errors in the synthetic NDVI images, but the MOD13Q1 product still has an important advantage over the MCD43A4, as the spatial resolution of MOD13Q1 is the highest among MODIS products. In addition, performing topographic correction on the Landsat and MODIS images may improve the results of STARFM outputs, especially in mountainous regions.

The spectral difference between pixels in base pair Landsat and MODIS images, regarded as a measure of the homogeneity of the MODIS pixel, is one of the three factors that decide the weight function in STARFM [21]. Thus, spectral differences caused by other reasons, such as different sensors and spectral response function, will lower the accuracy of synthetic images or increase the differences of the mean values between synthetic and MODIS images [33]. Therefore, the normalization of Landsat reflectance images could minimize these kinds of errors to obtain a consistent synthetic time series with MODIS.

In the study of Bhandari et al. [33] and Schmidt et al. [35], the number of base Landsat images used to generate the long-term synthetic time series were much larger than the number used in this study, with 38 scenes in five years and 90 scenes in 7.5 years, respectively. We only found 32 Landsat scenes in 12 years in this study area and, due to the unbalanced dates of observation and many scenes from the non-growing period, only six of them were finally used to generate the synthetic time series in our study. When using different base Landsat images with STARFM, Walker et al. [32] reported minor differences between the synthetic images on the same date. Therefore, despite the small amount of base Landsat images in our study, the accuracy of our synthetic time series were not significantly decreased. In addition, all of the six base Landsat images used in this study were obtained within the same season (between 21 August and 14 October), ensuring the same accuracy level of synthetic images at the same date in each year, which is important for the SMK test.

4.2. Potentials of the Trends Derived from Synthetic Time Series

At a spatial resolution of 250 m, NDVI trends derived from MODIS time series only captured the main structures of vegetation changes, while the fine-scale patterns of negative and positive trends were neutralized or even eliminated (Figures 6B,D and 9B,D). In contrast, at a spatial resolution of 30 m, synthetic time series could capture finer patterns of positive and negative trends at the sub-MODIS pixel scale in Figure 6A,C and Figure 9A,C. These results were supported by the study of Stellmes et al. [8], who examined the effects of spatial resolution on NDVI trends by spatially degraded Landsat images. Positive and negative slopes were neutralized in the trend map derived from MODIS time series, due to the large pixel size, while these two different slopes could be separated in the synthetic time series, thus reducing the proportion of pixels with a moderate slope value and enlarging the ranges (Figures 8 and 11). Therefore, synthetic time series captured vegetation change information that could not be captured by MODIS time series, making the distributions of NDVI slopes statistically different (p < 0.001) between the two time series. Consequently, STARFM-generated synthetic NDVI time series have the potential to enhance the ability to capture finer scale patterns of vegetation changes.

According to the statistical significance test (SMK and MK), the trends calculated from synthetic dense NDVI time series and NDVImax time series were different. Specifically, 68% of the vegetated areas were characterized by positive trend (41% being significant and 27% non-significant) in the synthetic dense time series, while the figure for MODIS dense time series was even higher, at 87% (56% being significant and 31% non-significant). By contrast, in the NDVImax time series, areas characterized by positive and negative trends were roughly at the same level, with 53% experiencing a positive trend (9% being significant and 43% non-significant) in synthetic series and 52% for MODIS series (5% being significant and 48% non-significant). The magnitudes of both positive and negative trends in NDVImax time series are higher than those in dense time series. These indicate that different vegetation trends were signified with different time series data around the Datong Coalfield. For example, an area may experience a positive trend in NDVI time series for the spring, autumn and winter seasons, but a negative one for the summer season, which may result in an increasing overall trend by examining dense NDVI time series with SMK test and TS slope for this area, but a decreasing trend by NDVImax time series with MK test and TS slope. Lasanta and Vicente-Serrano [9] reported different trend directions between August and March time series in northeast Spain, which supports the above assumption. This may explain the different trend directions between the dense NDVI time series and NDVImax time series described above. Consequently, it is worthy to perform trend analysis for NDVI time series of particular seasons or months, which may not be reflected in dense time series or NDVImax time series. STARFM also provides a feasible way to do such a job at a spatial resolution of 30 m, which is especially valuable in regions where the number of available Landsat scenes is very limited. However, field works need to be conducted to validate the trends obtained from STAFRM-generated time series data and to decide which time series is the most appropriate one on quantifying mining effects on vegetation.

STARFM usually generates poorer synthetic images in heterogeneous areas than in homogeneous areas [21], thus a lower consistency level or even totally different directions of the trends may be obtained between synthetic and MODIS time series in heterogeneous areas. The potential use of generating high spatial and high temporal resolution synthetic NDVI time series is to better quantify the effects of mining activities on vegetation by comparing vegetation trends between mined and non-mined areas. Therefore, the time series data generated in those areas should be accurate enough (i.e., consistent with MODIS time series) to obtain reliable results. Areas characterized by significantly different trends between synthetic and MODIS time series should not be selected to conduct such comparison work. An example would be the central part of the mined area, where the decreasing NDVI trend was obvious in the synthetic dense time series in Figure 6A,C, while an opposite NDVI trend was seen in the MODIS dense time series in Figure 6B,D. However, the decreasing NDVImax trend of this area was consistent in both the synthetic and MODIS time series (Figure 9), indicating that the synthetic NDVImax time series data were able to obtain consistent results in this area. This may be because errors in the synthetic images were reduced by calculating the annual maximum NDVI images.

Synthetic images generated by the same base Landsat scene were highly inter-correlated (Table 1, Figure 5), which was also reported by Schmidt et al. [35]. Therefore, the scarcity of Landsat scenes in the study area added the strong serial dependence of the synthetic time series, as they were produced by STARFM using base pair images from the adjacent year. Thus, the statistical methods that are robust against serial dependence used in this study were critical to obtain reliable NDVI trends from the synthetic dense time series.

5. Conclusions

Extensive coal mining activities pose a threat to the local environmental conditions of the Loess Plateau, where the ecosystems are very fragile. To accurately quantify the vegetation changes caused by coal mining activates, this study proposed a method to generate high spatial and high temporal time series NDVI data using the Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM). The main conclusions can be summarized as follows:

(1)

Using Landsat and MODIS NDVI images directly as input for the STARFM algorithm generated more accurate synthetic NDVI images than using band reflectance images as input for the STARFM algorithm around the Datong Coalfield.

(2)

It is feasible to use one base pair image to generate synthetic time series for a year and use base pair images from the adjacent years when no Landsat scenes are available in a certain year.

(3)

NDVI trends from synthetic dense time series (2000–2011) were found to be highly consistent with that from MODIS time series around the Datong Coalfield, while the consistency level was even higher for annual maximum NDVI (NDVImax) trends from synthetic and MODIS time series.

(4)

Despite the fact that the overall patterns of NDVI trends were maintained, trends from MODIS time series at a spatial resolution of 250 m neutralize fine-scale positive and negative vegetation changes. Synthetic time series generated via STARFM are likely to be more capable in capturing such fine-scale vegetation changes.

(5)

It is suggested that when using the synthetic time series to quantify the effects of mining activities on vegetation, the test areas should be selected with caution, as the NDVI trends obtained from synthetic high spatial resolution time series and MODIS time series may be significantly different in some heterogeneous areas.

Acknowledgments

This research was funded by the Research and Innovation Program for College and University Graduate Students in Jiangsu Province (CXZZ13_0938), the National Natural Science Foundation of China (Grant No. 51174287) and the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD). The NASA Goddard Space Flight Center and US Geological Survey (USGS) are thanked for sharing MODIS and Landsat data. We also wish to thank the Landsat Ecosystem Disturbance Adaptive Processing System for providing the STARFM and LEDAPS tools. We are grateful to the reviewers for their constructive and detailed comments and suggestions.

Conflict of Interest

The authors declare no conflict of interest.

References

  1. Comprehensive Management Planning Framework for Loess Plateau Area (2010–2030). Available online: http://www.sdpc.gov.cn/zcfb/zcfbtz/2010tz/W020110117531609590135.pdf (accessed on 29 August 2013).
  2. Laflen, J.M.; Tian, J.; Huang, C. Soil Erosion and Dryland Farming; CRC Press: Boca Raton, FL, USA, 2000; pp. 1–54. [Google Scholar]
  3. Lei, S.; Bian, Z.; Zhang, R.; Li, L. Study on the change rule of groundwater level and its impacts on vegetation at arid mining area. J. Coal Sci. Eng 2007, 13, 179–182. [Google Scholar]
  4. Zhang, J.Y.; Guo, L.W.; Gong, J.L.; Zhu, L.Q.; Gong, X.M. Analysis on impacting the soil environment of the underground coal mining. Prog. Min. Sci. Saf. Tech 2007, 8, 1966–1969. [Google Scholar]
  5. Zhang, X. Research on Temporal and Spatial Variation for Land Cover in Mining Disturbed Zone (In Chinese). M.Sc.,. Taiyuan University of Technology, Taiyuan, China, May 2010. [Google Scholar]
  6. Al-Hamdan, M.; Cruise, J.; Rickman, D.; Quattrochi, D. Effects of spatial and spectral resolutions on fractal dimensions in forested landscapes. Remote Sens 2010, 2, 611–640. [Google Scholar]
  7. Röder, A.; Udelhoven, T.; Hill, J.; del Barrio, G.; Tsiourlis, G. Trend analysis of Landsat-TM and -ETM+ imagery to monitor grazing impact in a rangeland ecosystem in Northern Greece. Remote Sens. Environ 2008, 112, 2863–2875. [Google Scholar]
  8. Stellmes, M.; Udelhoven, T.; Röder, A.; Sonnenschein, R.; Hill, J. Dryland observation at local and regional scale—Comparison of Landsat TM/ETM+ and NOAA AVHRR time series. Remote Sens. Environ 2010, 114, 2111–2125. [Google Scholar]
  9. Lasanta, T.; Vicente-Serrano, S.M. Complex land cover change processes in semiarid Mediterranean regions: An approach using Landsat images in northeast Spain. Remote Sens. Environ 2012, 124, 1–14. [Google Scholar]
  10. Parent, M.B.; Verbyla, D. The browning of Alaska’s boreal forest. Remote Sens 2010, 2, 2729–2747. [Google Scholar]
  11. Matejicek, L.; Kopackova, V. Changes in Croplands as a result of large scale mining and the associated impact on food security studied using time-series landsat images. Remote Sens 2010, 2, 1463–1480. [Google Scholar]
  12. Trejo, I.; Dirzo, R. Deforestation of seasonally dry tropical forest: a national and local analysis in Mexico. Biol. Conserv 2000, 94, 133–142. [Google Scholar]
  13. Stockli, R.; Vidale, P.L. European plant phenology and climate as seen in a 20-year AVHRR land-surface parameter dataset. Int. J. Remote Sens 2004, 25, 3303–3330. [Google Scholar]
  14. Slayback, D.A.; Pinzon, J.E.; Los, S.O.; Tucker, C.J. Northern hemisphere photosynthetic trends 1982–1999. Global Change Biol 2003, 9, 1–15. [Google Scholar]
  15. Julien, Y.; Sobrino, J.A. Global land surface phenology trends from GIMMS database. Int. J. Remote Sens 2009, 30, 3495–3513. [Google Scholar]
  16. Maignan, F.; Breon, F.M.; Bacour, C.; Demarty, J.; Poirson, A. Interannual vegetation phenology estimates from global AVHRR measurements—Comparison with in situ data and applications. Remote Sens. Environ 2008, 112, 496–505. [Google Scholar]
  17. Bai, Z.G.; Dent, D.L.; Olsson, L.; Schaepman, M.E. Proxy global assessment of land degradation. Soil Use Manag 2008, 24, 223–234. [Google Scholar]
  18. Fensholt, R.; Rasmussen, K.; Kaspersen, P.; Huber, S.; Horion, S.; Swinnen, E. Assessing land degradation/recovery in the African Sahel from long-term earth observation based primary productivity and precipitation relationships. Remote Sens 2013, 5, 664–686. [Google Scholar]
  19. Mao, J.; Shi, X.; Thornton, P.; Hoffman, F.; Zhu, Z.; Myneni, R. Global latitudinal-asymmetric vegetation growth trends and their driving mechanisms: 1982–2009. Remote Sens 2013, 5, 1484–1497. [Google Scholar]
  20. de Beurs, K.M.; Henebry, G.M. Trend Analysis of the Pathfinder AVHRR Land (PAL) NDVI—Data for the Deserts of Central Asia. IEEE Geosci. Remote Sens. Lett 2004, 1, 282–286. [Google Scholar]
  21. Gao, F.; Masek, J.; Schwaller, M.; Hall, F. On the blending of the Landsat and MODIS surface reflectance: Predicting daily Landsat surface reflectance. IEEE Trans. Geosci. Remote Sens 2006, 44, 2207–2218. [Google Scholar]
  22. Bechtel, B.; Zakšek, K.; Hoshyaripour, G. Downscaling land surface temperature in an urban area: A case study for Hamburg, Germany. Remote Sens 2012, 4, 3184–3200. [Google Scholar]
  23. Gao, F.; Kustas, W.; Anderson, M. A data mining approach for sharpening thermal satellite imagery over land. Remote Sens 2012, 4, 3287–3319. [Google Scholar]
  24. Wu, P.; Shen, H.; Ai, T.; Liu, Y. Land-surface temperature retrieval at high spatial and temporal resolutions based on multi-sensor fusion. Int. J. Digit. Earth 2013. [Google Scholar] [CrossRef]
  25. Zurita-Milla, R.; Clevers, J.G.P.W.; Schaepman, M.E. Unmixing-based Landsat TM and MERIS FR data fusion. IEEE Geosci. Remote Sens. Lett 2008, 5, 453–457. [Google Scholar]
  26. Wu, M.Q.; Niu, Z.; Wang, C.Y.; Wu, C.Y.; Wang, L. Use of MODIS and Landsat time series data to generate high-resolution temporal synthetic Landsat data using a spatial and temporal reflectance fusion model. J. Appl. Remote Sens. 2012, 6, 063507-1–063507-13. [Google Scholar]
  27. Acerbi-Junior, F.W.; Clevers, J.; Schaepman, M.E. The assessment of multi-sensor image fusion using wavelet transforms for mapping the Brazilian Savanna. Int. J. Appl. Earth Obs 2006, 8, 278–288. [Google Scholar]
  28. Hilker, T.; Wulder, M.A.; Coops, N.C.; Linke, J.; McDermid, G.; Masek, J.G.; Gao, F.; White, J.C. A new data fusion model for high spatial- and temporal-resolution mapping of forest disturbance based on Landsat and MODIS. Remote Sens. Environ 2009, 113, 1613–1627. [Google Scholar]
  29. Wulder, M.A.; Masek, J.G.; Cohen, W.B.; Loveland, T.R.; Woodcock, C.E. Opening the archive: How free data has enabled the science and monitoring promise of Landsat. Remote Sens. Environ 2012, 122, 2–10. [Google Scholar]
  30. Hilker, T.; Wulder, M.A.; Coops, N.C.; Seitz, N.; White, J.C.; Gao, F.; Masek, J.G.; Stenhouse, G. Generation of dense time series synthetic Landsat data through data blending with MODIS using a spatial and temporal adaptive reflectance fusion model. Remote Sens. Environ 2009, 113, 1988–1999. [Google Scholar]
  31. Singh, D. Evaluation of long-term NDVI time series derived from Landsat data through blending with MODIS data. Atmosfera 2012, 25, 43–63. [Google Scholar]
  32. Walker, J.J.; de Beurs, K.M.; Wynne, R.H.; Gao, F. Evaluation of Landsat and MODIS data fusion products for analysis of dryland forest phenology. Remote Sens. Environ 2012, 117, 381–393. [Google Scholar]
  33. Bhandari, S.; Phinn, S.; Gill, T. Preparing Landsat image time series (LITS) for monitoring changes in vegetation phenology in queensland, Australia. Remote Sens 2012, 4, 1856–1886. [Google Scholar]
  34. Roy, D.P.; Ju, J.; Lewis, P.; Schaaf, C.; Gao, F.; Hansen, M.; Lindquist, E. Multi-temporal MODIS–Landsat data fusion for relative radiometric normalization, gap filling, and prediction of Landsat data. Remote Sens. Environ 2008, 112, 3112–3130. [Google Scholar]
  35. Schmidt, M.; Udelhoven, T.; Gill, T.; Röder, A. Long term data fusion for a dense time series analysis with MODIS and Landsat imagery in an Australian Savanna. J. Appl. Remote Sens 2012, 6, 063512-1–063512-18. [Google Scholar]
  36. Meng, J.; Du, X.; Wu, B. Generation of high spatial and temporal resolution NDVI and its application in crop biomass estimation. Int. J. Digit. Earth 2011, 6, 203–218. [Google Scholar]
  37. Shen, H.; Wu, P.; Liu, Y.; Ai, T.; Wang, Y.; Liu, X. A spatial and temporal reflectance fusion model considering sensor observation differences. Int. J. Remote Sens 2013, 34, 4367–4383. [Google Scholar]
  38. Zhu, X.; Chen, J.; Gao, F.; Chen, X.; Masek, J.G. An enhanced spatial and temporal adaptive reflectance fusion model for complex heterogeneous regions. Remote Sens. Environ 2010, 114, 2610–2623. [Google Scholar]
  39. Jönsson, P.; Eklundh, L. TIMESAT—A program for analyzing time-series of satellite sensor data. Comput. Geosci 2004, 30, 833–845. [Google Scholar]
  40. Zhao, H.; Wang, Y.J. Research on the factors affecting the classification accuracy of ETM remote sensing image land cover/use (In Chinese). Remote Sens. Tech. Appl 2012, 27, 600–608. [Google Scholar]
  41. Fensholt, R.; Proud, S.R. Evaluation of earth observation based global long term vegetation trends—Comparing GIMMS and MODIS global NDVI time series. Remote Sens. Environ 2012, 119, 131–147. [Google Scholar]
  42. Xiao, X.M.; Braswell, B.; Zhang, Q.Y.; Boles, S.; Frolking, S.; Moore, B. Sensitivity of vegetation indices to atmospheric aerosols: Continental-scale observations in Northern Asia. Remote Sens. Environ 2003, 84, 385–392. [Google Scholar]
  43. Fensholt, R.; Nielsen, T.T.; Stisen, S. Evaluation of AVHRR PAL and GIMMS 10-day composite NDVI time series products using SPOT-4 vegetation data for the African continent. Int. J. Remote Sens 2006, 27, 2719–2733. [Google Scholar]
  44. Masek, J.G.; Vermote, E.F.; Saleous, N.E.; Wolfe, R.; Hall, F.G.; Huemmrich, K.F.; Gao, F.; Kutler, J.; Lim, T.K. A Landsat surface reflectance dataset for North America, 1990–2000. IEEE Geosci. Remote Sens. Lett 2006, 3, 68–72. [Google Scholar]
  45. Kaufman, Y.J.; Wald, A.E.; Remer, L.A.; Gao, B.-C.; Li, R.-R.; Flynn, L. The MODIS 2.1-μm channel-correlation with visible reflectance for use in remote sensing of aerosol. IEEE Trans. Geosci. Remote Sens 1997, 35, 1286–1298. [Google Scholar]
  46. Sen, P.K. Estimates of the regression coefficient based on Kendall’s Tau. J. Am. Stat. Assoc 1968, 63, 1379–1389. [Google Scholar]
  47. Theil, H. A rank-invariant method of linear and polynomial regression analysis. I, II and III. Nederl. Akad. Wetensch. Proc 1950, 53. [Google Scholar]
  48. Rousseeuw, P.J.; Leroy, A.M. Robust Regression and Outlier Detection; John Wiley & Sons, Inc: Hoboken, NJ, USA, 2003. [Google Scholar]
  49. Hirsch, R.M.; Slack, J.R. A nonparametric trend test for seasonal data with serial dependence. Water Resour. Res 1984, 20, 727–732. [Google Scholar]
  50. Vanbelle, G.; Hughes, J.P. Nonparametric-tests for trend in water-quality. Water Resour. Res 1984, 20, 127–136. [Google Scholar]
  51. de Beurs, K.M.; Henebry, G.M. A statistical framework for the analysis of long image time series. Int. J. Remote Sens 2005, 26, 1551–1573. [Google Scholar]
  52. Hirsch, R.M.; Slack, J.R.; Smith, R.A. Techniques of Trend Analysis for Monthly Water Quality Data. Water Resour. Res 1982, 18, 107–121. [Google Scholar]
  53. de Jong, R.; de Bruin, S.; de Wit, A.; Schaepman, M.E.; Dent, D.L. Analysis of monotonic greening and browning trends from global NDVI time-series. Remote Sens. Environ 2011, 115, 692–702. [Google Scholar]
  54. Fensholt, R.; Langanke, T.; Rasmussen, K.; Reenberg, A.; Prince, S.D.; Tucker, C.; Scholes, R.J.; Le, Q.B.; Bondeau, A.; Eastman, R.; et al. Greenness in semi-arid areas across the globe 1981–2007—An Earth Observing Satellite based analysis of trends and drivers. Remote Sens. Environ 2012, 121, 144–158. [Google Scholar]
  55. Yin, H.; Udelhoven, T.; Fensholt, R.; Pflugmacher, D.; Hostert, P. How Normalized Difference Vegetation Index (NDVI) trendsfrom Advanced Very High Resolution Radiometer (AVHRR) and Système Probatoire d’Observation de la Terre VEGETATION (SPOT VGT) time series differ in agricultural areas: An inner mongolian case study. Remote Sens 2012, 4, 3364–3389. [Google Scholar]
  56. Mann, H.B. Nonparametric tests against trend. Econometrica 1945, 13, 245–259. [Google Scholar]
  57. Kendall, M.G. Rank Correlation Methods, 4th ed; Charles Griffin: London, UK, 1975. [Google Scholar]
  58. Wang, X.L.L.; Swail, V.R. Changes of extreme wave heights in Northern Hemisphere oceans and related atmospheric circulation regimes. J. Climate 2001, 14, 2204–2221. [Google Scholar]
  59. Mirkin, B. Eleven ways to look at the chi-squared coefficient for contingency tables. Am. Stat 2001, 55, 111–120. [Google Scholar]
  60. Bhandari, S.; Phinn, S.; Gill, T. Assessing viewing and illumination geometry effects on the MODIS vegetation index (MOD13Q1) time series: implications for monitoring phenology and disturbances in forest communities in Queensland, Australia. Int. J. Remote Sens 2011, 32, 7513–7538. [Google Scholar]
Remotesensing 05 04255f1 1024
Figure 1. Location of the study area. The boundary of the Loess Plateau was provided by the National Development and Reform Commission. The lower right figure shows the Landsat image using an RGB (red-green-blue) composite of near infrared (NIR), red and green. The upper right figure shows the digital elevation model (DEM) obtained from the International Scientific Data Service Platform ( http://datamirror.csdb.cn/).
Figure 1. Location of the study area. The boundary of the Loess Plateau was provided by the National Development and Reform Commission. The lower right figure shows the Landsat image using an RGB (red-green-blue) composite of near infrared (NIR), red and green. The upper right figure shows the digital elevation model (DEM) obtained from the International Scientific Data Service Platform ( http://datamirror.csdb.cn/).
Remotesensing 05 04255f1 1024
Remotesensing 05 04255f2 1024
Figure 2. Land cover classification map of the study area (modified/adopted from Zhao et al. [40]).
Figure 2. Land cover classification map of the study area (modified/adopted from Zhao et al. [40]).
Remotesensing 05 04255f2 1024
Remotesensing 05 04255f3 1024
Figure 3. Two schemes for generation of synthetic Normalized Difference Vegetation Index (NDVI) images using Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM). t0 represents the day on which base pair Landsat and Moderate Resolution Imaging Spectroradiometer (MODIS) images were obtained; tn represents the day on which synthetic image was produced.
Figure 3. Two schemes for generation of synthetic Normalized Difference Vegetation Index (NDVI) images using Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM). t0 represents the day on which base pair Landsat and Moderate Resolution Imaging Spectroradiometer (MODIS) images were obtained; tn represents the day on which synthetic image was produced.
Remotesensing 05 04255f3 1024
Remotesensing 05 04255f4 1024
Figure 4. Per-pixel comparison between Landsat NDVI of 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000 and the corresponding synthetic NDVI generated from STARFM with base pair images of (A) NIR and red band reflectance of 13 October 2000 (scheme 1), (B) NDVI data of 13 October 2000 (scheme 2) and (C) NDVI data of 11 October 2002 (scheme 2, adjacent year). The times of the synthetic NDVI were the MODIS composite periods.
Figure 4. Per-pixel comparison between Landsat NDVI of 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000 and the corresponding synthetic NDVI generated from STARFM with base pair images of (A) NIR and red band reflectance of 13 October 2000 (scheme 1), (B) NDVI data of 13 October 2000 (scheme 2) and (C) NDVI data of 11 October 2002 (scheme 2, adjacent year). The times of the synthetic NDVI were the MODIS composite periods.
Remotesensing 05 04255f4 1024
Remotesensing 05 04255f5 1024
Figure 5. (A) Mean (black solid lines) and standard deviation (gray areas) values of the images in synthetic and MODIS dense NDVI time series (2000–2011). (B) Mean (black circles and red dots) and standard deviation (vertical bars) values of the images in synthetic and MODIS NDVImax time series (2000–2011).
Figure 5. (A) Mean (black solid lines) and standard deviation (gray areas) values of the images in synthetic and MODIS dense NDVI time series (2000–2011). (B) Mean (black circles and red dots) and standard deviation (vertical bars) values of the images in synthetic and MODIS NDVImax time series (2000–2011).
Remotesensing 05 04255f5 1024
Remotesensing 05 04255f6 1024
Figure 6. Spatial distribution of NDVI slopes (in NDVI units per year) derived from (A) synthetic and (B) MODIS dense NDVI time series (2000–2011) and the corresponding significance of the trends according to the seasonal Mann-Kendall test for (C) synthetic and (D) MODIS time series (2000–2011).
Figure 6. Spatial distribution of NDVI slopes (in NDVI units per year) derived from (A) synthetic and (B) MODIS dense NDVI time series (2000–2011) and the corresponding significance of the trends according to the seasonal Mann-Kendall test for (C) synthetic and (D) MODIS time series (2000–2011).
Remotesensing 05 04255f6 1024
Remotesensing 05 04255f7 1024
Figure 7. Scatterplots of the significance (Z scores) of NDVI trends versus the magnitude of NDVI trends (NDVI slopes) derived from (A) synthetic and (B) MODIS dense NDVI time series (2000–2011). Pixels beyond the two vertical dashed lines of Z = −1.96 and Z = 1.96 are characterized by a significant trend. N is the total number of vegetated pixels in the study area.
Figure 7. Scatterplots of the significance (Z scores) of NDVI trends versus the magnitude of NDVI trends (NDVI slopes) derived from (A) synthetic and (B) MODIS dense NDVI time series (2000–2011). Pixels beyond the two vertical dashed lines of Z = −1.96 and Z = 1.96 are characterized by a significant trend. N is the total number of vegetated pixels in the study area.
Remotesensing 05 04255f7 1024
Remotesensing 05 04255f8 1024
Figure 8. Histograms of the frequency of the NDVI slopes derived from (A) synthetic and (B) MODIS dense time series (2000–2011).
Figure 8. Histograms of the frequency of the NDVI slopes derived from (A) synthetic and (B) MODIS dense time series (2000–2011).
Remotesensing 05 04255f8 1024
Remotesensing 05 04255f9 1024
Figure 9. Spatial distribution of NDVImax slopes (in NDVI units per year) derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011) and the corresponding significance of the trends according to the Mann-Kendall test for (C) synthetic and (D) MODIS NDVImax time series.
Figure 9. Spatial distribution of NDVImax slopes (in NDVI units per year) derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011) and the corresponding significance of the trends according to the Mann-Kendall test for (C) synthetic and (D) MODIS NDVImax time series.
Remotesensing 05 04255f9 1024
Remotesensing 05 04255f10 1024
Figure 10. Scatterplots of the significance (Z scores) of NDVImax trends versus the magnitude of NDVI trends (NDVImax slopes) derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011). Pixels beyond the two vertical dashed lines of Z = −1.96 and Z = 1.96 are characterized by a significant trend. N is the total number of vegetated pixels in the study area.
Figure 10. Scatterplots of the significance (Z scores) of NDVImax trends versus the magnitude of NDVI trends (NDVImax slopes) derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011). Pixels beyond the two vertical dashed lines of Z = −1.96 and Z = 1.96 are characterized by a significant trend. N is the total number of vegetated pixels in the study area.
Remotesensing 05 04255f10 1024
Remotesensing 05 04255f11 1024
Figure 11. Histograms of the frequency of the NDVImax slopes derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011).
Figure 11. Histograms of the frequency of the NDVImax slopes derived from (A) synthetic and (B) MODIS NDVImax time series (2000–2011).
Remotesensing 05 04255f11 1024
Table
Table 1. Base pair images used in the generation of synthetic Landsat NDVI time series. MOD13Q1 products with a composite period closest to the observation time of Landsat images were selected as the base MODIS images. No available Landsat images were found in this study area for years highlighted (bold font), and base pair images from the adjacent year were used as inputs for STARFM.
Table 1. Base pair images used in the generation of synthetic Landsat NDVI time series. MOD13Q1 products with a composite period closest to the observation time of Landsat images were selected as the base MODIS images. No available Landsat images were found in this study area for years highlighted (bold font), and base pair images from the adjacent year were used as inputs for STARFM.
YearDate of Base LandsatComposite Period of Base MODIS
200010/13/200009/29∼10/14/2000
200108/21/200108/13∼08/28/2001
200209/25/200209/22∼10/07/2002
200309/25/200209/22∼10/07/2002
200410/14/200610/08∼10/23/2006
200510/14/200610/08∼10/23/2006
200610/14/200610/08∼10/23/2006
200709/15/200709/06∼09/21/2007
200809/15/200709/06∼09/21/2007
200910/06/200909/14∼09/29/2009
201010/06/200909/14∼09/29/2009
201110/06/200909/14∼09/29/2009
Table
Table 2. R2 and mean absolute difference among different vegetation types between Landsat NDVI images of 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000 and the corresponding synthetic NDVI images, with base pair images of NIR and red band reflectance of 13 October 2000 (scheme 1), NDVI data of 13 October 2000 (scheme 2) and NDVI data of 11 October 2002 (scheme 2, adjacent year).
Table 2. R2 and mean absolute difference among different vegetation types between Landsat NDVI images of 6 May 2000, 1 July 2000, 11 September 2000 and 21 October 2000 and the corresponding synthetic NDVI images, with base pair images of NIR and red band reflectance of 13 October 2000 (scheme 1), NDVI data of 13 October 2000 (scheme 2) and NDVI data of 11 October 2002 (scheme 2, adjacent year).
Base InputsVegetation TypesR2Mean Absolute Difference
05/0607/0109/1110/2305/0607/0109/1110/23
NIR and red band reflectance of 13 October 2000 (scheme 1)Non-irrigated0.350.450.260.410.0200.0530.0490.024
Irrigated0.340.670.490.460.0370.0530.0570.040
Shrub0.440.610.610.530.0410.0680.0490.041
Broadleaf0.540.430.590.420.0630.0850.0520.062
Grassland0.320.570.470.360.0260.0420.0450.031
NDVI data of 13 October 2000 (scheme 2)Non-irrigated0.470.520.370.510.0120.0440.0360.014
Irrigated0.540.720.550.770.0200.0480.0500.015
Shrub0.730.750.770.650.0270.0530.0380.026
Broadleaf0.810.550.760.700.0390.0660.0410.029
Grassland0.610.700.560.600.0130.0330.0380.013
NDVI data of 11 October 2000 (scheme 2, adjacent year)Non-irrigated0.460.510.310.600.0150.0450.0370.016
Irrigated0.410.730.520.570.0250.0510.0540.022
Shrub0.650.690.700.780.0220.0600.0390.017
Broadleaf0.680.380.570.670.0540.0730.0560.040
Grassland0.480.700.550.530.0130.0340.0410.012
Table
Table 3. Cross-tabulation considering the significance of NDVI trends derived from the synthetic and MODIS dense time series (2000–2011). Chi-square test of the cross-tabulation showed a significant (p < 0.001) relationship between the two trend maps with a contingency coefficient of 0.56 out of 0.86 (values in km2).
Table 3. Cross-tabulation considering the significance of NDVI trends derived from the synthetic and MODIS dense time series (2000–2011). Chi-square test of the cross-tabulation showed a significant (p < 0.001) relationship between the two trend maps with a contingency coefficient of 0.56 out of 0.86 (values in km2).
MODIS
Negative (p < 0.05)Negative (p > 0.05)Positive (p < 0.05)Positive (p > 0.05)Percentages
Synthetic
Negative (p < 0.05)258.350.4308.1120.515%
Negative (p > 0.05)46.559.5517.4191.117%
Positive (p < 0.05)53.341.31,556.9288.641%
Positive (p > 0.05)48.555.9274.0880.827%
Percentages9%4%56%31%100%
Table
Table 4. Cross-tabulation considering the significance of NDVImax trends derived from the synthetic and MODIS time series (2000–2011). The chi-square test of the cross-tabulation showed a significant (p < 0.001) relationship between the two trend maps with a contingency coefficient of 0.57 out of 0.86 (values in km2).
Table 4. Cross-tabulation considering the significance of NDVImax trends derived from the synthetic and MODIS time series (2000–2011). The chi-square test of the cross-tabulation showed a significant (p < 0.001) relationship between the two trend maps with a contingency coefficient of 0.57 out of 0.86 (values in km2).
MODIS
Negative (p < 0.05)Negative (p > 0.05)Positive (p < 0.05)Positive (p > 0.05)Percentages
Synthetic
Negative (p < 0.05)262.8155.87.196.112%
Negative (p > 0.05)90.5862.324.5580.836%
Positive (p < 0.05)9.839.5118.1214.69%
Positive (p > 0.05)55.2534.272.41218.343%
Percentages10%37%5%48%100%

Share and Cite

MDPI and ACS Style

Tian, F.; Wang, Y.; Fensholt, R.; Wang, K.; Zhang, L.; Huang, Y. Mapping and Evaluation of NDVI Trends from Synthetic Time Series Obtained by Blending Landsat and MODIS Data around a Coalfield on the Loess Plateau. Remote Sens. 2013, 5, 4255-4279. https://doi.org/10.3390/rs5094255

AMA Style

Tian F, Wang Y, Fensholt R, Wang K, Zhang L, Huang Y. Mapping and Evaluation of NDVI Trends from Synthetic Time Series Obtained by Blending Landsat and MODIS Data around a Coalfield on the Loess Plateau. Remote Sensing. 2013; 5(9):4255-4279. https://doi.org/10.3390/rs5094255

Chicago/Turabian Style

Tian, Feng, Yunjia Wang, Rasmus Fensholt, Kun Wang, Li Zhang, and Yi Huang. 2013. "Mapping and Evaluation of NDVI Trends from Synthetic Time Series Obtained by Blending Landsat and MODIS Data around a Coalfield on the Loess Plateau" Remote Sensing 5, no. 9: 4255-4279. https://doi.org/10.3390/rs5094255

APA Style

Tian, F., Wang, Y., Fensholt, R., Wang, K., Zhang, L., & Huang, Y. (2013). Mapping and Evaluation of NDVI Trends from Synthetic Time Series Obtained by Blending Landsat and MODIS Data around a Coalfield on the Loess Plateau. Remote Sensing, 5(9), 4255-4279. https://doi.org/10.3390/rs5094255

Article Metrics

Back to TopTop