Next Article in Journal
Analysis of Thermal Anomalies in Volcanic Areas Using Multiscale and Multitemporal Monitoring: Vulcano Island Test Case
Next Article in Special Issue
Estimation of Carbon Fluxes from Eddy Covariance Data and Satellite-Derived Vegetation Indices in a Karst Grassland (Podgorski Kras, Slovenia)
Previous Article in Journal
Urban Tomographic Imaging Using Polarimetric SAR Data
Previous Article in Special Issue
Generation of High Resolution Vegetation Productivity from a Downscaling Method
 
 
Correction published on 27 May 2019, see Remote Sens. 2019, 11(10), 1246.
Retraction published on 27 June 2019, see Remote Sens. 2019, 11(13), 1522.
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

RETRACTED: Estimation of Vegetation Productivity Using a Landsat 8 Time Series in a Heavily Urbanized Area, Central China

1
Research Center of Forest Remote Sensing and Information Engineering, Central South University of Forestry & Technology, Changsha 410004, China
2
Key Laboratory of Forestry Remote Sensing Based Big Data and Ecological Security for Hunan Province, Changsha 410004, China
3
Key Laboratory of State Forestry Administration on Forest Resources Management and Monitoring in Southern Area, Changsha 410004, China
*
Authors to whom correspondence should be addressed.
Remote Sens. 2019, 11(2), 133; https://doi.org/10.3390/rs11020133
Submission received: 15 October 2018 / Revised: 2 January 2019 / Accepted: 6 January 2019 / Published: 11 January 2019
(This article belongs to the Special Issue Remote Sensing of Primary Productivity)

Abstract

:
Estimating the net primary production (NPP) of vegetation is essential for eco-environment conservation and carbon cycle research. Remote sensing techniques, combined with algorithm models, have been proven to be promising methods for NPP estimation. High-precision and real-time NPP monitoring in heterogeneous areas requires high spatio-temporal resolution remote sensing data, which are not easy to acquire by single remote sensors, especially in cloudy weather. This study proposes to fuse images of different sensors to provide high spatio-temporal resolution data for NPP estimation in cloud-prone areas. Firstly, the time series Normalized Difference Vegetation Index (NDVI) with a spatial resolution of 30 m and a temporal resolution of 16 days, are obtained by the enhanced spatial and temporal adaptive reflectance fusion model (ESTARFM). Then, the time series NDVI data, combined with meteorological data are input into an improved Carnegie–Ames–Stanford Approach (CASA) model for NPP estimation. This method is validated by a case study of a heavily urbanized area, in the middle reaches of the Yangtze River in China. The results indicate that the NPP estimated by the fused NDVI data has more detailed spatial information than by using the MODIS data. The results show a strong correlation between the actual Landsat8 NDVI and the fused NDVI images, which means that the accuracy of synthetic NDVI images (a 16 day interval and a 30 m resolution) is reliable, and it can provide superior inputs for accurate estimations of a NPP time series. The correlation coefficient (R) and root mean square error between the NPP, based on the fused NDVI and the measured NPP, are 0.66 and 14.280 g C/(m2·yr), respectively, indicating a good consistency. The small discrepancy is caused by the uncertainties of fused NDVI, measurement errors, conversion errors, and other factors in the CASA model. In this study, we achieved NPP with high spatial and temporal resolutions, which can provide higher accuracies of NPP data for analyzing the carbon cycling heavily urbanized areas, compared with similar studies using mono-temporal NPP data. The spatio-temporal fusion technique is an effective way of generating high spatio-temporal resolution images from different sensors, thereby providing enough data for NPP monitoring in urbanized areas.

1. Introduction

Vegetation is an important factor for the regulation of CO2 exchange between the atmosphere and the land surface [1,2,3]. Urban vegetation is an important sink for carbon cycling in urban ecosystems, which can also reduce fossil fuel usage through the processes of transpiration, shading, and the blocking of winds [4,5,6]. Therefore, measurements of urban vegetation carbon storage can help with the understanding of the relationship between urban vegetation and global carbon, accounting for CO2 emissions and improved urban planning and management [7,8,9].
Net primary productivity (NPP) can be defined as the total amount of carbon accumulated as vegetation biomass [10,11], which is a fundamental function of terrestrial ecosystems. NPP is an important indicator of ecosystem function, and a key element in the carbon cycle, since it shows the carbon storage capacity of a plant community, as well as the vegetation’s response to climate and land cover changes [12,13,14]. Recently, land cover changes, such as the extension of urban areas, have brought about many environmental problems, including the reduction of biodiversity, urban heat island effect enhancement, environment pollution aggravation, and the disorder of carbon circulation. Therefore, timely and accurate estimation and monitoring of NPP is of great significance for regional ecosystem assessments and carbon cycle research.
Remote sensing technology can provide information on vegetation dynamics, which is useful for carbon budget estimation. By being able to offer data of the same area at different times, remote sensing technology has become an important method of NPP estimation at different spatial and temporal scales [15,16,17]. Based on satellite sensors, a variety of models has been developed to estimate NPP, such as the Carnegie–Ames–Stanford Approach (CASA) [18] and other models [19,20,21,22,23,24,25]. For example, the CASA model uses a light-use efficiency (LUE) factor, which is the efficiency of the conversion of light energy into dry materials by vegetation [18,26]. The CASA model is frequently employed to estimate NPP with different types of remote sensing data [27,28,29].
The NPP of large-scale areas (e.g., global, continental, and national) has been estimated at a spatial resolution of 1.1 km by the Advanced Very High Resolution Radiometer data, and at 500 m by using Moderate Resolution Imaging Spectrometer (MODIS) data [30,31,32]. MODIS data operated by the National Aeronautics and Space Administration of the United States (NASA) are widely used, due to their moderate temporal and spatial resolutions, as well as free access to data [33,34]. The spatial resolution for NPP modeling will influence the accuracy of NPP estimation. Regions with temperate climates have a continuous gradation in land cover types [35,36,37]. Heavily urbanized regions have different spatial patterns, coupled with highly fragmented vegetation cover. Thus, remote sensing data at high spatial/temporal resolutions are important for accurate NPP estimation.
With higher spatial resolution, Landsat TM/ETM+/OLI can provide high accuracy NPP estimations for regions with highly fragmented vegetation. However, the long revisit intervals of sensors, and rainy weather and cloud coverage, limit the number of usable cloud-free images. Integrating images from different sensors is a feasible and affordable way of achieving a time series of images with high spatio-temporal resolutions. With the development of remote sensing technology, many spatial and temporal fusion models have been proposed [38,39]. The spatial and temporal adaptive reflectance fusion model (STARFM) proposed by Gao et al. (2006) [40], performs well in preserving spatial and temporal details. Based on this, several improved models have been developed, such as the spatial temporal adaptive algorithm for mapping reflectance change (STAARCH) [41], the enhanced STARFM (ESTARFM) [42], the operational STARFM data fusion framework, the spatio-temporal integrated temperature fusion model (STITFM) [43], and the robust adaptive spatial and temporal fusion model (RASTARFM) [44]. There are also some other spatial and temporal satellite image fusion models [45,46]. Among these fusion models, the most widely used STARFM model can capture temporal changes from coarse-resolution images, and still maintain the spatial details of the fine-resolution images [44]. However, STARFM does not consider the complexity of land cover changes in heterogeneous landscapes, resulting in the prediction of images with poor accuracy. The ESTARFM fusion model has demonstrated its potential in improving the accuracy of phenological change predictions in heterogeneous regions [42,44]. Thus, we used ESTARFM in this study to achieve a time series Landsat-like NDVI.
The object of this study is to develop a method for estimating NPP with high spatio–temporal resolution, using a fused NDVI time series and the CASA model. We choose the Changsha–Zhuzhou–Xiangtan area of Hunan Province, China, a humid and heavily urbanized region, as a study area. The MOD17A3H and field data are used to validate the accuracy of the estimated NPP.

2. Study Area and Data

2.1. Study Area

In this study, a heavily urbanized area in China was selected. The study area, the core area of the Changsha–Zhuzhou–Xiangtan area, is located south of the Yangtze River and north of Nanling Mountain, between the latitudes 27°36′N and 28°33′N, and between the longitudes 112°36′E and 113°16′E, in Hunan Province, China (Figure 1). It covers the urban areas of Changsha, Zhuzhou, and Xiangtan cities, as well as Changsha, Wangcheng, Zhuzhou, and Xiangtan counties; a total of 15 districts (counties) with an area of about 8640 km2. The variation of altitude in the study area ranges roughly from −7 m to 809 m. The landscape features mountains, hills, and inter-mountainous plains. The percentage of mountainous, hilly, and plain areas are about 14.7%, 28.2%, and 32.1%, respectively. The region has a subtropical monsoon climate with four distinct seasons. The average annual temperature is approximately 16–18 °C, and the annual precipitation is about 1814 mm. The vegetation types in the study area are mainly forest, grassland, and crops. The forest is distributed in the northeast, southeast, and southwest of the study area. The distribution of the grassland is mainly along the river bench. Crops mainly scatter around the lakes and rivers for the sake of irrigation.

2.2. Data and Processing

In this study, we used MOD13Q1 and Landsat-8 OLI to generate fused time series Landsat-like NDVI data. The Digital Elevation Model (DEM) data was used to obtain spatialized meteorological data, which, together with the fused NDVI time series and the Land use/Land cover data, were employed to generate the time series Landsat-like data. MOD17A3H and field data were utilized to validate the estimated NPP.

2.2.1. MODIS13Q1 Data

The MODIS NDVI was derived from the MODIS13Q1 products, with a spatial resolution of 250 m in the Sinusoidal projection [47,48]. The MODIS products data (23 MODIS13Q1 images) between January 2015 and December 2015, were acquired from the United States Geological Survey (USGS). After removing the invalid values by using the pixel reliability images, all of the MODIS NDVI time series were transformed into the UTM (WGS84) projection, zone 49 (North), the same as Landsat 8 OLI. Then, the pre-processed MODIS NDVI was smoothed by the Savizky–Golay (S-G) filter.

2.2.2. Landsat 8 OLI Data

The Landsat 8 OLI consists of eight multispectral bands (30 m), one panchromatic band (15 m), and two thermal bands (100 m). We assembled clear (0–5% cloud cover) Landsat 8 OLI images covering the study area in 2015 (including two adjacent overlapping Landsat footprints, path/row: 123/40 and 123/41) for the experiment (https://www.usgs.gov/) [49]. Radiometric calibration was performed to convert the digital number (DN) value to surface spectral reflectance, and atmospheric correction was conducted by using the FLAASH model of ENVI, version 5.1 [49]. The NDVI was computed by band 4 (red) and band 5 (NIR) using the surface spectral reflectance images [50,51].

2.2.3. DEM and Meteorological Data

A DEM with a spatial resolution of 30 m was used for the study (http://www.gscloud.cn/) [52]. The meteorological data, including the monthly mean temperature, monthly mean precipitation, and the monthly total solar radiation of the meteorological station (a total of 46 stations) in 2015, in and around the study area, were acquired from the China Meteorological Data Sharing Service System (http://cdc.nmic.cn) [52]. In order to improve the accuracy of the meteorological data, the multiple linear regression method was adopted to spatialize the meteorological elements [52]. First, we extracted the latitude, longitude, grade, and the slope grid layer of the study area at the pixel scale from the DEM. The grid layers, together with the DEM and the latitude and longitude information of the meteorological station, were then employed to obtain the spatialized meteorological data by a multiple linear regression method. The mean annual relative error between the observation data from the meteorological station and the spatialized meteorological data was less than 0.05, which meets the requirements of this study.

2.2.4. Land Use/Land Cover Data

Land use/land cover data were acquired from the National Earth System Science Data Sharing Infrastructure (http://www.geodata.cn) [52]. It provides a land cover remote sensing survey and monitoring database of China at a scale of 1:100,000. The dataset was generated by the visual interpretation of Landsat TM/ETM remote sensing images [53]. This product has six land use types, which can be further divided into 25 subclasses
In order to further analyze the NPP of several major land use types in the Changsha–Zhuzhou–Xiangtan area, we reclassified the land use types to five groups: forest, grassland, farmland, built-up (building pixels with vegetation), and other land (other land pixels with vegetation) (Figure 2).

2.2.5. The MOD17A3H NPP Product

The MOD17A3H annual NPP product in 2015 of the study area was selected as the validation data. The MOD17A3H Version 6 product provides information about the annual (yearly) NPP at a 500 m pixel resolution. MOD17A3H updates the Biome Property Look Up Tables of MOD17A3, and it uses the meteorological data of the updated Global Modeling and Assimilation Office [54]. MOD17A3H has a high accuracy, and it passed the test conducted by eddy covariance flux towers distributed in different climate zones and biozones [55,56]. The MODIS Reprojection Tool (MRT) was used to transform the MOD17A3H tiles into a UTM projection. Pixels with NPP values out of the range of the triple standard deviation, or with ineffective values, were omitted. The MOD17A3H data of the study area is shown in Figure 3.

2.2.6. Field Sampling Data

The aboveground biomass data used in this study were collected in October, 2015. Based on the investigated vegetation information, i.e., the vegetation growth and distribution of the study area, we selected 92 sites with geographical coordinates confirmed by Global Position System (GPS). Stratified sampling was used for sampling, and field sampling types include forest (pine, evergreen magnolia, and fir), grassland (sedge and clover), and farmland (cotton, paddy, and ramie). Each sample site had an area of 30 m × 30 m, representing a typical vegetation community with homogeneous vegetation and land cover types. In each site, five subsample squares (2 m × 2 m) were set up for vegetation collection. To obtain the actual aboveground biomass of the sample sites, all of the aboveground plants in the five plots (2 m × 2 m) of each sample site were harvested and dried. The dry weight of the biomass was then converted to the weight of carbon by using the carbon content, which was considered to be 45% in this study [57]. Measuring NPP for large areas is very difficult. However, the vegetation productivity and biomass are strongly correlated with NPP. Therefore, we applied the relationship between the aboveground biomass and NPP to validate the simulated NPP.

3. Methods

In this paper, we developed a method for estimating the NPP with high spatial–temporal resolution using a fused NDVI time series and an improved CASA model. The proposed method has three main steps, and a detailed flowchart is shown in Figure 4. First, the ESTARFM model [42] was employed to generate a synthetic NDVI with high spatial resolution and temporal frequency from Landsat and MODIS images. Second, we used an improved CASA model [58] to obtain NPP with high temporal and spatial resolutions, based on the acquired synthetic NDVI, as well as the meteorological and land use/cover data. Finally, we utilized the measured biomass to validate the synthetic NDVI-derived NPP data, based on the relationship between the aboveground biomass and NPP.

3.1. MODIS–NDVI Time-Series Filtering

The MODIS–NDVI data were generated from the MOD13Q1 data, using the maximum-value composite method, which reduces the noise caused by cloud and aerosol effects. Other noises and inaccurate phenologies were eliminated by an S-G filter, using TIMESAT software [59]. This filter could clearly describe minor changes in the study area, despite the complex crop types and broken plots. The S-G filtering was performed with a locally adapted moving window, which utilizes a polynomial least-squares regression to fit the time-series data as follows:
Y j = i = m i = m C i Y j + 1 N
where Y and Y* represent the original and fitted NDVIs of the time series, respectively, j is the running index of the original ordinate data table, Ci represents the coefficient of the ith NDVI value of the filter, N is the sliding array width, which is equal to the size of the smoothing window (2m + 1), and m is the half-width of the smoothing window [59]. In this study, the sizes of the moving window, the adaptive strength, and the number of envelope iterations were set as 5, 2, and 2, respectively, by the trial and error method. The effects of the S-G algorithm are shown in Figure 5, in which sudden abnormal drops in the NDVI time series have been adjusted exactly.

3.2. The ESTARFM Model for Synthetic NDVI

The ESTARFM model employs two pairs of fine- and coarse-resolution images obtained on the same date, and a coarse resolution image obtained on the prediction day [42]. The two fine-resolution images were used to search for pixels that are similar to the central pixel from the two coarse-resolution images in a moving window, and the pair of fine- and coarse-resolution images with better search results was then employed to synthesize the Landsat-like image [42]. During the fusion, the algorithm centers on a predicted pixel to set a domain window with a certain size. Then, convolution computation is carried out for all of the pixels in the window, using a weighting function W. Finally, the predicted value of the central pixel is obtained, using Equation (2):
L ( x w / 2 , y w / 2 , t p ) = L ( x w / 2 , y w / 2 , t 0 ) + i = 1 N W i × V i × ( M ( x w / 2 , y w / 2 , t p ) M ( x w / 2 , y w / 2 , t 0 )
where N is the number of similar pixels, w is the size of the moving window, (xi, yi) is the position of the ith similar pixel, and Wi is the weight of the ith similar pixel, which is determined by the spatial distance, time interval, and the spectral interval [42].
We selected the Landsat8 NDVI data from DOY164 and DOY212, and the MODIS NDVI data from DOY161 and DOY209 (the time nearest to the acquisition date of Landsat8 data) as the input base pair images. Then, we used the ESTARFM algorithm to predict time series Landsat-like images (16 day intervals and 30 m resolution). In order to examine the quality of the predicted NDVI images, actual Landsat8 OLI images were used as the reference to evaluate the fused images. The information of the input and the validating images of this study are listed in Table 1.

3.3. The NPP Estimation Model: An Improved CASA Model

The CASA model is robust in estimating NPP spatial and temporal patterns from remote sensing data. In this study, an improved CASA model was used to estimate NPP (Equation (3)) [58]. As illustrated in this model, weather variables and land cover data are two main input factors. The NPP (g C/(m2·month)) of pixel x in month t is the product of the absorbed photosynthetically active radiation (APAR) (MJ/(m2·month)) and the light use efficiency factor ε (g C/MJ).
N P P ( x , t ) = A P A R ( x , t ) × ε ( x , t )
A P A R ( x , t ) = S O L ( x , t ) × F P A R ( x , t ) × 0.5
ε ( x , t ) = ε m a x × T ( x , t ) × W ( x , t )
where SOL (x, t) (MJ/(m2·month)) is the total solar radiation of pixel x in month t. FPAR represents the absorption ratio of APAR by the vegetation layer. The constant of 0.5 represents the ratio of the effective solar radiation used by the vegetation to the total solar radiation (Equation (4)) [15]. εmax is the maximal light use efficiency under ideal conditions. It is related to the vegetation fractions and types. T and W are the temperature stress coefficient and moisture stress coefficient, respectively (Equation (5)) [58,60].
There are two differences between this model and the conventional CASA model. First, this model has simplified the calculation of the water stress coefficient. The spatial distribution of the water stress coefficient as calculated by Equations (6)–(8) is more consistent with the actual situation [58,60]:
W ( x , t ) = 0.5 + 0.5 × E ( x , t ) / E p ( x , t )
E ( x , t ) = { P ( x , t ) × R n ( x , t ) × [ ( P ( x , t ) ) 2 + ( R n ( x , t ) ) 2 + P ( x , t ) × R n ( x , t ) ] } / { [ P ( x , t ) + R n ( x , t ) ] × [ ( R n ( x , t ) ) 2 × ( R n ( x , t ) ) 2 ] }
E p ( x , t ) = [ E ( x , t ) + E p 0 ( x , t ) ] / 2
where E(x, t), P(x, t), and Rn(x, t) represent the actual evapotranspiration (ET), precipitation, and surface net radiation of pixel x in month t, respectively. Ep(x, t) and Ep0(x, t) are the potential ET and local potential ET.
Second, the potential maximum light utilization (εmax) is estimated by different types of vegetation, based on the simulation results of the BIOME–BGC model [61] and the measured NPP in the study area, which is more consistent with the actual situation [62].
The improved CASA model calculates the NPP by month, so we synthesized the 23 predicted NDVI into 12 months by the maximum value composite algorithm. Note that the predicted NDVI on DOY145 was used twice. In this study, NDVImax and NDVImin were calculated by using the fused time series and land use/cover data. The εmax in the traditional CASA model is 0.389 g C/MJ, which is always amended according to the specific vegetation types in the study area. The εmax of the forest and grassland, estimated based on the simulation results of the BIOME–BGCmodel [61] and the actual situation in China by Zhu et al. [63] and Pei et al. [64], are closer to the εmax derived from the measured NPP, compared with the conventional CASA model. According to Zhu et al. [58], farmland has the same or similar εmax as grassland. However, fertilization and irrigation could largely improve farmland productivity, especially in China, where farmers apply a lot of fertilizer and irrigation. Thus, the εmax values for farmland could be much higher than that of grassland. Crop yield is a part of NPP in the growing season in a certain region, and there is an effective yield conversion relationship between them [65]. Therefore, in this study, we used the NPP derived from crop yield as the measured NPP, and then the εmax value of farmland (0.61 g C/MJ) was simulated according to the error function formula. The NPP derived from crop yield can be calculated by the following equation:
N P P = i = 1 N Y i × M R Y i ( 1 M C i ) × 0.45 H I i × 0.9 / i = 1 N h a r v e s t e d a r e a i
where N, Yi, and MRYi are the number of crops, the yield of crop i, and the mass per unit of yield for crop i. MCi and HIi are the harvest indices and harvest moisture content of crop i, respectively [65]. The potential variability within the species for each of the factors that is used to convert yield to NPP was not considered in this study.
The parameters of the improved CASA model were confirmed based on the predicted time series NDVI data (16 day intervals and 30 m resolution) and the land cover data (Table 2).

3.4. Validation of the Estimated NPP

To assess the performance of the improved CASA algorithm with the time series MODIS NDVI (250 m) and the fused Landsat-like NDVI (30 m), MOD17A3H and field sampling data were used to verify the NPP estimates. We used 92 field sampling data to assess the accuracy of NPP derived from fused Landsat-like NDVI. During the validation of the NPP derived from MODIS data, 85 sample points were randomly selected in the study area. We used two statistical indices, the coefficient of correlation (R), and the root mean square error (RMSE). The RMSE is calculated as follows:
R M S E = i = 1 n ( y i y ¯ i ) 2 / n
where yi is the measured NPP or MOD17A3H value, y ¯ i is the mean of the estimated values, and n is the number of samples used in this study [15].

4. Result and Analysis

4.1. Accuracy of the Predicted NDVI

The MODIS NDVI, fused Landsat 8 NDVI, and actual Landsat 8 NDVI of the study area on DOY257 are shown in Figure 6. The fused NDVI image using ESTARFM had a much higher spatial resolution, and more details than the corresponding MODIS image. The fused NDVI and actual Landsat8 NDVI images were visually similar, even at the junctions of complex land covers, such as along the river bench, near the lake shore, and in some sparse vegetation areas. However, it was reported that the fused NDVI images using STARFM have always had over-estimations in areas near rivers or lakes, because of the abrupt surface changes caused by urban flooding [66,67]. Therefore, fused Landsat-like NDVI derived from ESATRFM was more reliable than that of SATRFM. Nevertheless, cloud contamination may cause bias between the actual Landsat8 NDVI and the fused NDVI images predicted by ESATRFM. Although the MODIS NDVI employed in this paper had a 16-day temporal resolution, it could not fully remove the noise in regions with cloudy and rainy weather [68]. In general, the fused NDVI (30 m) images not only contained high-frequency temporal information for MODIS–NDVI, but also the high spatial information of Landsat NDVI.
Figure 7 shows the scatterplots between the actual Landsat8 NDVI and the fused NDVI of DOY 257 of the pixel values. Most of the points are distributed around the 1:1 diagonal. The parameters of these two images are shown in Table 3. The results show that the fused image derived from ESTARFM has a high accuracy, and that the mean NDVI value and standard deviation of the actual image and the fused image are similar. However, the absolute errors in the STARFM predictions of NDVI tend to be the largest in the tree cover type, followed by shrubs, which have been demonstrated by Xie et al. [69]. Generally, correlations between the actual Landsat8 NDVI and the fused NDVI images are good. Therefore, the accuracies of the fused NDVI images (16 day intervals and 30 m resolution) generated by the ESTARFM algorithm were more reliable for NPP estimations than that of STARFM.

4.2. NPP Estimation Results and Validation

The NPP estimated from MODIS NDVI (Figure 8a) had a low spatial resolution and obscure expression in small areas, and it lacked detailed information in space. The NPP estimated from the fused NDVI (16 day intervals and 30 m resolution) had a higher spatial resolution (Figure 8b) and detailed spatial information, showing differences between small objects. The differences between the river, agricultural land, construction land, and forest were obvious, and the boundaries among them were clear.
The mean NPP of the study area in 2015 based on MODIS NDVI was 260.88 g C/(m2·yr), and the NPP based on the fused NDVI was 272.15 g C/(m2·yr), which was nearly 4.3% higher than the former. This small difference lies with the pixels in some sparse vegetation areas, for which the fused NDVI has a higher value than the MODIS NDVI [65]. Most pixel values of the NPP based on MODIS NDVI range between 200–600 g C/(m2·yr), where it is difficult to show the spatial heterogeneity of the land cover types. However, the pixel values of the NPP based on the fused NDVI data obey the normal distribution, with a peak value approximating 300 g C/(m2·yr), which is similar to that of a real distribution.
The results show that when the pixel vegetation coverage was less than 20%, the NPP value based on fused NDVI (30 m) was obviously larger than that based on the MODIS NDVI (Figure 9a). This is because the MODIS NDVI had a lower spatial resolution than the fused NDVI, so that the pixels with low vegetation coverage may be regarded as mixed pixels with non-vegetation coverage. The fused NDVI had a higher spatial resolution; so it had a smaller probability of forming mixed pixels with non-vegetation coverage. When the pixel vegetation coverage was larger than 80%, the NPP values based on these two NDVI were in good agreement (Figure 9b). The higher the pixel vegetation coverage was, the higher the probability of it being treated as a pure pixel. If both the MODIS NDVI and the fusion NDVI were considered as pure pixels, their estimated NPP values would be the same.

4.3. Accuracy Assessment of NPP Based on the Improved CASA Model

We first used the MOD17A3H data to validate the NPP values estimated by the CASA model using the MODIS data (250 m). The NPP derived from MODIS was resampled to 500 m, to be consistent with the MOD17A3H. The results are shown in Figure 10. The correlation coefficient and RMSE between the estimated NPP, based on MODIS NDVI and MOD17A3H, were nearly 0.87 and 15.427 g C/(m2·yr), indicating that the improved CASA model has a high level of accuracy.
In order to assess the accuracy of the NPP based on the fused NDVI data using an improved CASA model, we calculated the correlation between the field observation data and the estimated values in the study area. In this study, we used the relationship between the biomass and NPP to verify the accuracy of the simulated NPP, since the actual observation data was biomass. The correlation coefficient and the RMSE between the simulated NPP and measured NPP were 0.81 and 14.280 g C/(m2·yr), respectively (Figure 11). Therefore, the simulated and measured NPP had a significant linear relationship. The results also indicated that the estimated NPP based on fused NDVI using the improved CASA model was precise enough for further analysis.

4.4. Distribution and Variation of NPP using the Fused NDVI

The overall regional distribution of NPP in the study area indicates a high NPP in the northeast, southeast, and southwest, and a low NPP in the center of the Changsha–Zhuzhou–Xiangtan area, where a large number of built-up areas are distributed (Figure 8b). Figure 12 illustrates the spatial distribution of NPP in 2015 by month. NPP for April–September exhibited a large spatial heterogeneity, which is of great significance for detecting vegetation conditions in the study area. In April, the vegetation NPP was generally below 26 g C/(m2·yr), and in some parts, it was 22.1 g C/(m2·yr). During May, nearly half of the vegetation NPP values reached 55–60 g C/(m2·yr). The vegetation NPP values continued to increase from June to July, and they reached their peak between the late July and early August. In July, the highest value was 216 g C/(m2·yr) in the southeast region, but it was as low as 30 g C/(m2·yr) in the urban region. The vegetation NPP started decreasing during the last week of August. In the second week of September, the regional mean vegetation NPP was 30 g C/(m2·yr), and it decreased to 22.3 g C/(m2·yr) in October.
The estimated NPP time series describes more detailed spatial variations of NPP at a resolution of 30 m. The time series NPP retained not only the high-frequency temporal information from the MODIS NDVI time series, but also the high spatial information from Landsat. It was expected that such data can provide reliable data and scientific methods for the NPP changing analysis in heavily urbanized regions.
We then investigated the relationship between NPP variation and the land use/cover types, which is beneficial to evaluate the effects of land cover on the NPP quantity and structure. The mean NPP values of different land covers in the Changsha–Zhuzhou–Xiangtan area (Figure 2) in 2015 are shown in Figure 13. The average NPP of forest was 441.25 g C/(m2·yr), grassland was 289.61 g C/(m2·yr), farmland was 392.85 g C/(m2·yr), built-up with vegetation was 202.83 g C/(m2·yr), and other land with vegetation was 105.18 g C/(m2·yr). The spatial distribution of NPP was related to the land cover types. Mountains or hills in the northeast, southeast, and southwest of the study area were covered by extensive forests, resulting in a higher NPP. The central part of the study area was dominated by build-up zones, so that it had low NPP values (Figure 13a). The monthly NPP value of forest was relatively stable, with a continuous covering of broad-leaved forests during the whole year in most of the study area. Buildings and other land had nearly consistent monthly NPPs, because of the low vegetation coverage in these areas. The monthly NPP of grassland and farmland grew gradually from January to August, and then it reduced gradually. This trend was consistent with the plant growing trend in the study area (Figure 13b).

5. Discussion

NPP indicates the vegetation ecosystem service ability, and it reflects the variation of plant responses to climate change and human activities. However, land use/cover changes may directly alter the local ecological structure, resulting in NPP changes. A rapid decrease of carbon storage is observed in developing urban or coastal areas. Early studies mainly used MODIS data with a spatial resolution of 500 m or 250 m for urban NPP estimation. However, NPP derived from MODIS cannot meet precision NPP monitoring requirements in complex and heterogenous urban landscapes. In order to obtain high spatial resolution NPP estimations, some studies have resorted to moderate (Landsat TM/ETM+) or high spatial resolution remote sensing images [70,71,72]. However, these NPP data are either based on single-date images, or on a few scenes that could not show the changes of NPP in rapidly developing urban areas. In this study, NPP with high spatial and temporal resolutions in a heavily urbanized area was estimated based on the fused Landsat NDVI time series, which will be useful for monitoring and analyzing carbon cycling in urban areas.
High-precision and timely NPP monitoring in urbanized regions requires remote sensing data with high spatial–temporal resolutions, which are very difficult to acquire by single remote sensors, especially in urban areas with humid, cloudy climates. Integrating data from different sensors can achieve a time series with high temporal/spatial resolutions, providing appropriate data for regional vegetation dynamic monitoring and NPP estimation. The ESTARFM algorithm has been applied successfully in monitoring seasonal changes, and in capturing phenological information, as it can improve the accuracy of classifications for land cover. Phenologies extracted from the fused time series generated by ESATRFM using moderate resolution images are useful in paddy rice mapping [73,74]. The mean values of the fused NDVI obtained by the ESATRFM are almost identical to those of the MODIS NDVI, and the predicted NDVI time series keeps the high-frequency temporal information from the MODIS NDVI time series [75,76]. Therefore, the synthetic NDVI time series is reliable for estimating the NPP time series. In this study, we obtained the fused NDVI time series with a resolution of 16 days. Time series with higher temporal resolutions can be obtained if needed. The daily NDVI time series can improve the identification of critical growth phases of vegetation, but it would result in larger datasets and greater probabilities of cloud contamination.
We obtained the synthetic NDVI time series with a higher spatial resolution by ESTARFM; nevertheless, there were some uncertainties in the results of fused NDVI. A number of factors could affect the prediction accuracy in the synthetic images. The first factor is the selected base image pairs. Extra image pairs that are available during the same time period can improve the data fusion accuracy, especially on some key pair dates (e.g., the green-up and senescence stages for crops) [69]. We may be able to develop high-quality data fusion results for a long time series with a few key pairs of input. The second factor is the residual cloud contamination in the 16-day MODIS time series NDVI, which is quite common in tropical and sub-tropical areas [67]. Synthetic aperture radar (SAR) is a promising way to estimate the NPP time series, because of its well-timed image acquisitions, and the independence of its meteorological conditions. However, it is expensive to acquire a time series of SAR images covering large scale areas. The third factor is the land cover change in the study area. Although the ESATRFM can improve the accuracy of phenological change predictions in heterogeneous regions, it cannot accurately predict land cover changes [44]. The ESTARFM algorithm developed from SATRFM is useful only when Landsat data and MODIS data are both available at the same time. Some studies have applied ESTARFM and other fusion models to other satellite datasets [76]. Although the ESTARFM can derive fine spatial and temporal scales, it should be carefully handled for other datasets.
In this study, the improved CASA model proposed by Zhu et al. (2005) [58] was employed to estimate the NPP in the Changsha–Zhuzhou–Xiangtan area, and MOD17A3H and field measurements were used to validate the estimation results. The validation results show that the improved CASA model, combined with the fused NDVI time series, can achieve a high accuracy in the NPP estimation.
The value of εmax varies with different vegetation types, and its size is controversial, since it has a great influence on the estimation of NPP [77]. The improved CASA used the same εmax value for grassland and farmland. However, there are many studies pointing out that farmland has much higher εmax values than grassland [28,65,78]. Fertilization and irrigation could largely improve farmland productivity, especially in China. Thus, the εmax values for farmland could be much higher than that of grassland. In this study, we used the NPP derived from crop yield as the measured NPP to estimate the εmax value of farmland (0.61 g C/MJ). We assigned a εmax value above 0 to the built-up areas in this study. The εmax value of the same land cover type also varied with regions, such as built-up areas. The built-up pixels in the study area are usually mixed with vegetation, because green vegetation is always planted around or on top of many buildings, so that an εmax value that is larger than 0 was assigned. Zhang et al. [15] and Zhu et al. [58] also set the value of εmax to 0.542 and 0.389 for urban lands in southwest China and Inner Mongolia, respectively.
Considering the cost, coverage, and spatial resolution, Landsat 8 spectral images were selected, as they are usually available online free of charge, with a considerable image width and a long historical data record [79], which is useful for NPP time series estimations incorporating large areas and many years. The study shows that the proposed method is promising for NPP estimation, based on moderate spatial resolution remote sensing images. However, there are still some limitations for the proposed method in complex heterogeneous areas with heavy land cover changes. Pixel mixing is the first problem, because of the moderate spatial resolution and the complex heterogeneity of the study area. Higher resolution multi-spectral images (e.g., GF-1, GF-5) are a feasible method for the estimation of NPP with higher accuracy, but it also faces cost and image coverage problems [80]. In addition, mixed pixels of fused NDVI using spatial and temporal fusion models also need further investigation. The second factor is residual cloud contamination in the 16-day MODIS NDVI time series, which is quite common in tropical and sub-tropical areas. Optical images, combined with SAR images, are a promising method for NPP estimation, due to well-timed image acquisitions and the independence of SAR from meteorological conditions. The third factor is the method of NPP validation. In this study, we used field-measured biomass to validate the estimated NPP, which inevitably has some errors that are caused by measurement errors and conversion errors, resulting in some impact on verification results [15].
With accelerating urbanization, the land use/cover types of urban areas are experiencing dramatic changes, which also have a significant impact on the net primary productivity of vegetation. Therefore, a time series NPP of many years and of NPP spatio-temporal dynamics under land use/cover changes in urban areas need further analysis.

6. Conclusions

This study proposed a method to estimate the NPP in a heavily urbanized area by an improved CASA model, using the fused NDVI time series. This study has demonstrated the potential of using moderate spatial resolution images, combined with the CASA model, for NPP estimation of large areas. Despite the mixed pixel issues, issues during data fusion, and other potential issues affecting the estimation accuracy of NPP, the correlation coefficient and the RMSE between MODIS-derived NPP and MOD17A3H are 0.76 and 15.427 g C/(m2·yr), respectively, and the correlation coefficient and RMSE between NPP based on fused NDVI and field sample data are 0.66 and 14.280g C/(m2·yr), respectively. Therefore, the proposed method has the potential to provide an acceptable dataset of the spatial distribution of NPP analysis under land use/cover changes in urban areas.

Author Contributions

M.Z. designed and conducted the experiments and wrote the manuscript. H.L., G.W., H.S. provided suggestions for the analysis, discussion, and manuscript writing. Y.C. provided and processed the remote data.

Funding

This work was supported by the Twelfth Five-Year Plan Pioneering project of the High Technology Plan of the National Department of Science and Technology (2012AA102001), and the Scientific Research Fund of Hunan Provincial Education Department (17A225): Synergy Simulation of multi-resolution remote sensing data for city vegetation carbon mapping and uncertainty analysis.

Acknowledgments

The authors would like to thank the editors and anonymous reviewers for valuable comments, which are significant for improving this manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Turner, D.P.; Ritts, W.D.; Cohen, W.B.; Gower, S.T.; Running, S.W.; Zhao, M.S.; Costa, M.H.; Kirschbaum, A.A.; Ham, J.M.; Saleska, S.R.; et al. Evaluation of MODIS NPP and GPP products across multiple biomes. Remote Sens. Environ. 2006, 102, 282–292. [Google Scholar] [CrossRef]
  2. Crabtree, R.; Potter, C.; Mullen, R.; Sheldon, J.; Huang, S.L.; Harmsen, J.; Rodman, A.; Jean, C.A. Modeling and spatio-temporal analysis framework for monitoring environmental change using NPP as an ecosystem indicator. Remote Sens. Environ. 2009, 113, 1486–1496. [Google Scholar] [CrossRef]
  3. Bandaru, V.; West, T.O.; Ricciuto, D.M.; Izaurralde, R.C. Estimating crop net primary production using national inventory data and MODIS-derived parameters. ISPRS J. Photogramm. Remote Sens. 2013, 80, 61–71. [Google Scholar] [CrossRef]
  4. Wetherley, E.B.; Mcfadden, J.P.; Roberts, D.A. Megacity-scale analysis of urban vegetation temperatures. Remote Sens. Environ. 2018, 213, 18–33. [Google Scholar] [CrossRef]
  5. Myeong, S.; Nowak, D.J.; Duggin, M.J. A temporal analysis of urban forest carbon storage using remote sensing. Remote Sens. Environ. 2006, 101, 277–282. [Google Scholar] [CrossRef]
  6. Tigges, J.; Lakes, T.; Hostert, P. Urban vegetation classification: Benefits of multitemporal RapidEye satellite data. Remote Sens. Environ. 2013, 136, 66–75. [Google Scholar] [CrossRef]
  7. Lu, Y.; Coops, N.C.; Hermosilla, T. Estimating urban vegetation fraction across 25 cities in pan-Pacific using Landsat time series data. ISPRS J. Photogramm. Remote Sens. 2016, 126, 11–23. [Google Scholar] [CrossRef]
  8. Anchang, J.Y.; Ananga, E.O.; Pu, R. An efficient unsupervised index based approach for mapping urban vegetation from IKONOS imagery. Int. J. Appl. Earth Obs. Geoinf. 2016, 50, 211–220. [Google Scholar] [CrossRef]
  9. Rougier, S.; Puissant, A.; Stumpf, A.; Lachiche, N. Comparison of sampling strategies for object-based classification of urban vegetation from very high resolution satellite image. Int. J. Appl. Earth Obs. Geoinf. 2016, 51, 60–73. [Google Scholar] [CrossRef]
  10. Group on Earth Observations. Forest Carbon Tracking Portal. Available online: http://www.geo-fct.org/ (accessed on 23 April 2013).
  11. Fu, Y.; Lu, X.; Zhao, Y.; Zeng, X.; Xia, L. Assessment Impacts of Weather and Land Use/Land Cover (LULC) Change on Urban Vegetation Net Primary Productivity (NPP): A Case Study in Guangzhou, China. Remote Sens. 2013, 5, 4125–4144. [Google Scholar] [CrossRef]
  12. Piao, S.; Ciais, P.; Friedlingstein, P.; Peylin, P.; Reichstein, M.; Luyssaert, S.; Margolis, H.; Fang, J.; Barr, A.; Chen, A. Net carbon dioxide losses of northern ecosystems in response to autumn warming. Nature 2008, 451, 49–52. [Google Scholar] [CrossRef]
  13. Guan, X.; Shen, H.; Gan, W.; Yang, G.; Wang, L.; Li, X. A 33-year NPP monitoring study in southwest China by the fusion of multi-source remote sensing and station data. Remote Sens. 2017, 9, 1082. [Google Scholar] [CrossRef]
  14. Piao, S.; Fang, J.; He, J. Variations in vegetation net primary production in the Qinghai-Xizang Plateau, China, from 1982 to 1999. Clim. Chang. 2006, 74, 253–267. [Google Scholar] [CrossRef]
  15. Zhang, R.; Zhou, Y.; Luo, H.; Wang, F.; Wang, S. Estimation and Analysis of Spatiotemporal Dynamics of the Net Primary Productivity Integrating Efficiency Model with Process Model in Karst Area. Remote Sens. 2017, 9, 477. [Google Scholar] [CrossRef]
  16. Sellers, P.; Meeson, B.; Hall, F.; Asrar, G.; Murphy, R.; Schiffer, R.; Bretherton, F.; Dickinson, R.; Ellingson, R.; Field, C. Remote sensing of the land surface for studies of global change: Models—algorithms—experiments. Remote Sens. Environ. 1995, 51, 3–26. [Google Scholar] [CrossRef]
  17. Patenaude, G.; Hill, R.A.; Milne, R.; Gaveau, D.L.A.; Briggs, B.B.J.; Dawson, T.P. Quantifying forest aboveground carbon content using lidar remote sensing. Remote Sens. Environ. 2004, 93, 368–380. [Google Scholar] [CrossRef]
  18. Potter, C.S.; Randerson, J.T.; Field, C.B.; Matson, P.A.; Vitousek, P.M.; Mooney, H.A.; Klooster, S.A. Terrestrial ecosystem production: A process model-based on global satellite and surface data. Glob. Biogeochem. Cycles 1993, 7, 811–841. [Google Scholar] [CrossRef]
  19. Zhao, M.; Heinsch, F.A.; Nemani, R.R.; Running, S.W. Improvements of the MODIS terrestrial gross and net primary production global data set. Remote Sens. Environ. 2005, 95, 164–176. [Google Scholar] [CrossRef]
  20. Prince, S.D.; Goward, S.N. Global primary production: A remote sensing approach. J. Biogeogr. 1995, 22, 815–835. [Google Scholar] [CrossRef]
  21. Knorr, W.; Heimann, M. Impact of drought stress and other factors on seasonal land biosphere CO2 exchange studied through an atmospheric tracer transport model. Tellus 1995, 47, 471–489. [Google Scholar] [CrossRef]
  22. Sellers, P.J.; Randall, D.A.; Collatz, G.J. A revised land surface parameterization (SiB2) for atmospheric GCMs. I. Model formulation. J. Clim. 1996, 9, 676–705. [Google Scholar] [CrossRef]
  23. Ruimy, A.; Dedieu, G.; Saugier, B. TURC: A diagnostic model of continental gross primary productivity and net primary productivity. Glob. Biogeochem. Cycles 1996, 10, 269–286. [Google Scholar] [CrossRef]
  24. Liu, J.; Chen, J.M.; Cihlar, J. A Process-Based Boreal Ecosystem Productivity Simulator Using Remote Sensing Inputs. Remote Sens. Environ. 1997, 62, 158–175. [Google Scholar] [CrossRef]
  25. Parton, W.J.; Scurlock, J.M.O.; Ojıma, D.S. Observations and modeling of biomass and soil organic matter dynamics for the grassland biome worldwide. Glob. Biogeochem. Cycles 1993, 7, 785–809. [Google Scholar] [CrossRef]
  26. Rafique, R.; Zhao, F.; Rogier, D.J.; Zeng, N.; Ghassem, A. Global and Regional Variability and Change in Terrestrial Ecosystems Net Primary Production and NDVI: A Model-Data Comparison. Remote Sens. 2016, 8, 177. [Google Scholar] [CrossRef]
  27. Running, S.; Ramakrishna, R.; Nemani, F.; Heinsch, A.; Maosheng, Z.; Reeves, M.; Hashimoto, H. Acontinuous satellite-derived measure of global terrestrial primary production. BioScience 2004, 54, 547–560. [Google Scholar] [CrossRef]
  28. Donmez, C.; Berberoglu, S.; Curran, P.J. Modelling the current and future spatial distribution of NPP in a Mediterranean watershed. Int. J. Appl. Earth Obs. Geoinf. 2011, 13, 336–345. [Google Scholar] [CrossRef]
  29. Pachavo, G.; Murwira, A. Remote sensing net primary productivity (NPP) estimation with the aid of GIS modelled shortwave radiation (SWR) in a Southern African Savanna. Int. J. Appl. Earth Obs. Geoinf. 2014, 30, 217–226. [Google Scholar] [CrossRef]
  30. Bao, G.; Bao, Y.; Qin, Z.; Xin, X.; Bao, Y.; Bayarsaikan, S.; Zhou, Y.; Chuntai, B. Modeling net primary productivity of terrestrial ecosystems in the semi-arid climate of the Mongolian Plateau using LSWI-based CASA ecosystem model. Int. J. Appl. Earth Obs. Geoinf. 2016, 46, 84–93. [Google Scholar] [CrossRef]
  31. Hasenauer, H.; Petritsch, R.; Zhao, M.; Boisvenus, C.; Running, S.W. Reconciling satellite with ground data to estimate forest productivity at national scales. For. Ecol. Manag. 2012, 276, 196–208. [Google Scholar] [CrossRef]
  32. Neumann, M.; Zhao, M.; Kindermann, G.; Hasenauer, H. Comparing MODIS Net Primary Production Estimates with Terrestrial National Forest Inventory Data in Austria. Remote Sens. 2015, 7, 3878–3906. [Google Scholar] [CrossRef]
  33. Bala, G.; Joshi, J.; Chaturvedi, R.K.; Gangamani, H.V.; Hashimoto, H.; Nemani, R. Trends and variability of AHVRR-derived NPP in India. Remote Sens. 2013, 5, 810–829. [Google Scholar] [CrossRef]
  34. Zhu, L.; Southworth, J. Disentangling the relationships between net primary production and precipitation in southern africa savannas using satellite observations from 1982 to 2010. Remote Sens. 2013, 5, 3803–3825. [Google Scholar] [CrossRef]
  35. Ichii, K.; Kondo, M.; Okabe, Y.; Ueyama, M.; Kobayashi, H.; Lee, S.-J.; Saigusa, N.; Zhu, Z.; Myneni, R.B. Recent changes in terrestrial gross primary productivity in Asia from 1982 to 2011. Remote Sens. 2013, 5, 6043–6062. [Google Scholar] [CrossRef]
  36. Xiao, X.; Hollinger, D.; Aber, J.; Goltz, M.; Davidson, E.A.; Zhang, Q.; Moore, B., III. Satellite-based modeling of gross primary production in an evergreen needle leaf forest. Remote Sens. Environ. 2004, 89, 519–534. [Google Scholar] [CrossRef]
  37. Ardö, J.; Tagesson, T.; Jamali, S.; Khatir, A. MODIS EVI-based net primary production in the Sahel 2000–2014. Int. J. Appl. Earth Obs. Geoinf. 2018, 65, 35–45. [Google Scholar] [CrossRef]
  38. Chen, B.; Huang, B.; Xu, B. Comparison of Spatiotemporal Fusion Models: A Review. Remote Sens. 2015, 7, 1798–1835. [Google Scholar] [CrossRef]
  39. Zhang, H.; Huang, B.; Zhang, M.; Yu, L. A generalization of spatial and temporal fusion methods for remotely sensed surface parameters. Int. J. Remote Sens. 2015, 36, 4411–4445. [Google Scholar] [CrossRef]
  40. 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]
  41. Hiker, H.; 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 dada blending with MODIS using a spatial and temporal adaptive reflectance fusion model. Remote Sens. Environ. 2009, 113, 1988–1999. [Google Scholar] [CrossRef]
  42. Zhu, X.L.; Chen, J.; Gao, F.; Chen, X.H.; 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] [CrossRef]
  43. Wu, P.; Shen, H.; Zhang, L.; Göttsche, F.M. Integrated fusion of multi-scale polar-orbiting and geostationary satellite observations for the mapping of high spatial and temporal resolution land surface temperature. Remote Sens. Environ. 2015, 156, 169–181. [Google Scholar] [CrossRef]
  44. Zhao, Y.; Huang, B.; Song, H. A robust adaptive spatial and temporal image fusion model for complex land surface changes. Remote Sens. Environ. 2018, 208, 42–62. [Google Scholar] [CrossRef]
  45. Wu, M.Q.; Niu, Z.; Wang, C.; Wu, C.; 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. [Google Scholar] [CrossRef]
  46. Wang, Q.; Atkinson, P.M. Spatio-temporal fusion for daily Sentinel-2 images. Remote Sens. Environ. 2017, 204, 65–80. [Google Scholar] [CrossRef]
  47. Levin, N.; Heimowitz, A. Mapping spatial and temporal patterns of Mediterranean wildfires from MODIS. Remote Sens. Environ. 2012, 126, 12–26. [Google Scholar] [CrossRef]
  48. Liao, C.; Wang, J.; Pritchard, I.; Liu, J.; Shang, J. A Spatio-Temporal Data Fusion Model for Generating NDVI Time Series in Heterogeneous Regions. Remote Sens. 2017, 9, 1125. [Google Scholar] [CrossRef]
  49. Mansaray, L.; Huang, W.; Zhang, D.; Huang, J.; Li, J. Mapping Rice Fields in Urban Shanghai, Southeast China, Using Sentinel-1A and Landsat 8 Datasets. Remote Sens. 2017, 9, 257. [Google Scholar] [CrossRef]
  50. Hasituya; Chen, Z.; Wang, L.; Wu, W.; Jiang, Z.; Li, H. Monitoring Plastic-Mulched Farmland by Landsat-8 OLI Imagery Using Spectral and Textural Features. Remote Sens. 2016, 8, 353. [Google Scholar] [CrossRef]
  51. Karakizi, C.; Karantzalos, K.; Vakalopoulou, M.; Antoniou, G. Detailed land cover mapping from multitemporal landsat-8 data of different cloud cover. Remote Sens. 2018, 8, 1214. [Google Scholar] [CrossRef]
  52. Chen, X.; Zeng, Y. Spatial and temporal variability of the net primary production (NPP) and its relationship with climate factors in subtropical mountainous and hilly regions of China: A case study in Hunan province. Acta Geogr. Sin. 2016, 71, 35–48. [Google Scholar]
  53. Zhou, Y.; Xiao, X.; Qin, Y.; Dong, J.; Zhang, G.; Kou, W.; Wang, J.; Li, X. Mapping paddy rice planting area in rice-wetland coexistent areas through analysis of Landsat 8 OLI and MODIS images. Int. J. Appl. Earth Obs. 2016, 46, 1–12. [Google Scholar] [CrossRef] [PubMed]
  54. Cracknell, A.P.; Kanniah, K.D.; Tan, K.P.; Wang, L. Towards the development of a regional version of MOD17 for the determination of gross and net primary productivity of oil palm trees. Int. J. Remote Sens. 2015, 36, 262–289. [Google Scholar] [CrossRef]
  55. Malone, S.L.; Tulbure, M.G.; Perez-Luque, A.J.; Assal, T.; Bremer, L.; Drucker, D.P.; Hillis, V.; Varela, S.; Goulden, M. Drought resistance across California ecosystems: Evaluating changes in carbon dynamics using satellite imagery. Ecosphere 2016, 7, e01561. [Google Scholar] [CrossRef]
  56. Running, S.; Mu, Q.; Zhao, M. MOD17A3H MODIS/Terra Net Primary Production Yearly L4 Global 500 m SIN Grid V006 Data Set; NASA EOSDIS Land Processes DAAC: Washington, DC, USA, 2015.
  57. Johnson, W.C.; Sharpe, D.M. The ratio of total to merchantable forest biomass and its application to the global carbon budget. Can. J. For. Res. 1983, 13, 372–383. [Google Scholar] [CrossRef]
  58. Zhu, W.; Pan, Y.; Hu, H.; Li, J.; Gong, P. Estimating net primary productivity of terrestrial vegetation based on remote sensing: A case study in Inner Mongolia, China. J. Remote Sens. 2005, 9, 300–307. [Google Scholar]
  59. Jönsson, P.; Eklundh, L. A program for analyzingtime-series of satellite sensor data. Comput. Geosci. 2004, 30, 833–845. [Google Scholar] [CrossRef]
  60. Zhou, G.; Zhang, X. A natural vegetation NPP model. Acta Phytoecol. Sin. 1995, 19, 193–200. [Google Scholar]
  61. Running, S.W.; Thornton, P.E.; Nemani, R.; Glassy, J.M. Global Terrestrial Gross and Net Primary Productivity from the Earth Observing System. Methods in Ecosystem Science; Springer: New York, NY, USA, 2000; pp. 44–57. [Google Scholar]
  62. Zhu, W. Estimation of Net Primary Productivity of Chinese Terrestrial Vegetation Based on Remote Sensing and Its Relationship with Global Climate Change; Beijing Normal University: Beijing, China, 2005. [Google Scholar]
  63. Zhu, W.; Pan, Y.; Hao, H.; Yu, D.; Hu, H. Simulation of maximum light use efficiency for some typical vegetation types in China. Chin. Sci. Bull. 2006, 51, 457–463. [Google Scholar] [CrossRef]
  64. Pei, F.; Li, X.; Liu, X.; Wang, S.; He, Z. Assessing the differences in net primary productivity between pre- and post-urban land development in China. Agric. For. Meteorol. 2013, 171–172, 174–186. [Google Scholar] [CrossRef]
  65. Lobell, D.B.; Hicke, J.A.; Asner, G.P.; Field, C.B.; Tucker, C.J.; Los, S.O. Satellite estimates of productivity and light use efficiency in United States agriculture, 1982-98. Glob. Chang. Biol. 2002, 8, 722–735. [Google Scholar] [CrossRef]
  66. Zhang, F.; Zhu, X.L.; Liu, D.S. Blending MODIS and Landsat images for urban flood mapping. Int. J. Remote Sens. 2014, 35, 3237–3253. [Google Scholar] [CrossRef]
  67. Zhang, B.; Zhang, L.; Xie, D.; Yin, X.; Liu, C.; Liu, G. Application of Synthetic NDVI Time Series Blended from Landsat and MODIS Data for Grassland Biomass Estimation. Remote Sens. 2016, 8, 10. [Google Scholar] [CrossRef]
  68. Rembold, F.; Atzberger, C.; Savin, I.; Rojas, O. Using low resolution satellite imagery for yield prediction and yield anomaly detection. Remote Sens. 2013, 5, 1704–1733. [Google Scholar] [CrossRef]
  69. Xie, D.; Gao, F.; Sun, L.; Anderson, M. Improving Spatial-Temporal Data Fusion by Choosing Optimal Input Image Pairs. Remote Sens. 2018, 10, 1142. [Google Scholar] [CrossRef]
  70. Sun, R.; Chen, J.M.; Zhu, Q.; Zhou, Y.; Liu, J.; Li, J.; Liu, S.; Yan, G.; Tang, S. Spatial distribution of net primary productivity and evapotranspiration in Changbaishan Natural Reserve, China, using Landsat ETM+ data. Can. J. Remote Sens. 2004, 30, 731–742. [Google Scholar] [CrossRef]
  71. Wang, P.; Sun, R.; Hu, J.; Zhu, Q.; Zhou, Y.; Li, L.; Chen, J.M. Measurements and simulation of forest leaf area index and net primary productivity in Northern China. J. Environ. Manag. 2007, 85, 607–615. [Google Scholar] [CrossRef] [PubMed]
  72. Chen, J.M.; Chen, X.; Ju, W. Effects of vegetation heterogeneity and surface topography on spatial scaling of net primary productivity. Biogeosciences 2013, 10, 4879–4896. [Google Scholar] [CrossRef]
  73. Singha, M.; Wu, B.; Zhang, M. An Object-Based Paddy Rice Classification Using Multi-Spectral Data and Crop Phenology in Assam, Northeast India. Remote Sens. 2016, 8, 479. [Google Scholar] [CrossRef]
  74. Liu, M.; Liu, X.; Wu, L.; Zou, X.; Jiang, T.; Zhao, B. A Modified Spatiotemporal Fusion Algorithm Using Phenological Information for Predicting Reflectance of Paddy Rice in Southern China. Remote Sens. 2018, 10, 772. [Google Scholar] [CrossRef]
  75. Heimhuber, V.; Tulbure, M.G.; Broich, M. Addressing spatio-temporal resolution constraints in Landsat and MODIS-based image mapping of large-scale floodplain inundation dynamics. Remote Sens. Environ. 2018, 211, 307–320. [Google Scholar] [CrossRef]
  76. Gärtner, P.; Förster, M.; Kleinschmit, B. The benefit of synthetically generated RapidEye and Landsat 8 data fusion time series for riparian forest disturbance monitoring. Remote Sens. Environ. 2016, 177, 237–247. [Google Scholar] [CrossRef]
  77. Field, C.B.; Randerson, J.T.; Malmström, C.M. Global net primary production: Combining ecology and remote sensing. Remote Sens. Environ. 1995, 51, 74–88. [Google Scholar] [CrossRef]
  78. Bradford, J.B.; Hicke, J.A.; Lauenroth, W.K. The relative importance of light-use efficiency modifications from environmental conditions and cultivation for estimation of large-scale net primary productivity. Remote Sens. Environ. 2005, 96, 246–255. [Google Scholar] [CrossRef]
  79. Veh, G.; Korup, O.; Roessner, S.; Walz, A. Detecting Himalayan glacial lake outburst floods from Landsat time series. Remote Sens. Environ. 2018, 207, 84–97. [Google Scholar] [CrossRef]
  80. Zhou, D.; Zhao, S.; Zhang, L.; Liu, S. Remotely sensed assessment of urbanization effects on vegetation phenology in China’s 32 major cities. Remote Sens. Environ. 2016, 176, 272–281. [Google Scholar] [CrossRef]
Figure 1. Location of the Changsha–Zhuzhou–Xiangtan area.
Figure 1. Location of the Changsha–Zhuzhou–Xiangtan area.
Remotesensing 11 00133 g001
Figure 2. The reclassified land use/land cover map.
Figure 2. The reclassified land use/land cover map.
Remotesensing 11 00133 g002
Figure 3. MOD17A3H data of the study area.
Figure 3. MOD17A3H data of the study area.
Remotesensing 11 00133 g003
Figure 4. The flowchart of the proposed method for net primary production (NPP) estimation.
Figure 4. The flowchart of the proposed method for net primary production (NPP) estimation.
Remotesensing 11 00133 g004
Figure 5. A pixel of the Moderate Resolution Imaging Spectrometer–Normalized Difference Vegetation Index (MODIS–NDVI) time series after filtering by Savitzky-Golay (S-G).
Figure 5. A pixel of the Moderate Resolution Imaging Spectrometer–Normalized Difference Vegetation Index (MODIS–NDVI) time series after filtering by Savitzky-Golay (S-G).
Remotesensing 11 00133 g005
Figure 6. (a) The MODIS NDVI of DOY 257, (b) the fused NDVI generated by the ESTARFM model of DOY 257, and (c) the actual Landsat8 OLI NDVI of DOY 260.
Figure 6. (a) The MODIS NDVI of DOY 257, (b) the fused NDVI generated by the ESTARFM model of DOY 257, and (c) the actual Landsat8 OLI NDVI of DOY 260.
Remotesensing 11 00133 g006
Figure 7. (a) Scatterplots between the actual Landsat8 NDVI and the fused NDVI images of DOY257, and (b) scatterplots between the actual Landsat8 NDVI and MODIS NDVI.
Figure 7. (a) Scatterplots between the actual Landsat8 NDVI and the fused NDVI images of DOY257, and (b) scatterplots between the actual Landsat8 NDVI and MODIS NDVI.
Remotesensing 11 00133 g007
Figure 8. NPP estimations based on (a) MODIS NDVI data and (b) the fused NDVI data.
Figure 8. NPP estimations based on (a) MODIS NDVI data and (b) the fused NDVI data.
Remotesensing 11 00133 g008
Figure 9. Estimated NPP results under vegetation coverage (a) <20% and (b) >80%.
Figure 9. Estimated NPP results under vegetation coverage (a) <20% and (b) >80%.
Remotesensing 11 00133 g009
Figure 10. Comparison between the estimated NPP and MOD17A3H NPP.
Figure 10. Comparison between the estimated NPP and MOD17A3H NPP.
Remotesensing 11 00133 g010
Figure 11. Comparison between the estimated NPP and the measured NPP.
Figure 11. Comparison between the estimated NPP and the measured NPP.
Remotesensing 11 00133 g011
Figure 12. Monthly time series NPP of the Changsha–Zhuzhou–Xiangtan area in 2015, estimated based on the fused NDVI.
Figure 12. Monthly time series NPP of the Changsha–Zhuzhou–Xiangtan area in 2015, estimated based on the fused NDVI.
Remotesensing 11 00133 g012
Figure 13. NPP of different land cover types using the fused NDVI. (a) The annual gross NPP value and (b) the monthly NPP value.
Figure 13. NPP of different land cover types using the fused NDVI. (a) The annual gross NPP value and (b) the monthly NPP value.
Remotesensing 11 00133 g013
Table 1. Input images of the enhanced spatial and temporal adaptive reflectance fusion (ESTARFM) model for prediction and validation.
Table 1. Input images of the enhanced spatial and temporal adaptive reflectance fusion (ESTARFM) model for prediction and validation.
Input Landsat t1Input Landsat t2Input MODIS t1Input MODIS t2Input MODIS tkValidation Landsat tk
TestDOY164DOY212DOY161DOY209DOY257DOY260
Table 2. Land cover types and NPP model parameters in the study area.
Table 2. Land cover types and NPP model parameters in the study area.
Level 1Type CodeLevel 2NDVImaxNDVIminεmax (g C/MJ)
Forest11Evergreen coniferous forest0.7510.5620.39
12Evergreen broad-leaved forest0.7350.5290.97
13Deciduous coniferous forest0.7430.5580.51
14Deciduous broad-leaved forest0.7380.5520.66
15Coniferous and broad-leaved mixed forest0.7420.5760.49
16Shrub0.7530.6210.45
Grassland21Meadow grassland0.7120.3240.54
22Typical grassland0.7250.5360.54
26Shrub grassland0.7490.6350.54
Farmland31Paddy field0.7130.5340.61
32Irrigated land0.7320.5090.61
33Dry land0.7250.5280.61
Build-up41Urban construction land0.5320.2570.54
42Rural construction land0.6620.4230.54
Other land51Swamp0.6580.4120.54
53Inland water------
54River beach0.6350.2620.54
61Bare rock------
62Bare land------
Table 3. Accuracy of the MODIS NDVI (Figure 6a), the fused NDVI images (Figure 6b), and the actual Landsat8 NDVI (Figure 6c).
Table 3. Accuracy of the MODIS NDVI (Figure 6a), the fused NDVI images (Figure 6b), and the actual Landsat8 NDVI (Figure 6c).
DateDOY257 2015
Figure 6aFigure 6bFigure 6c
Mean value0.6950.6820.713
Standard deviation0.0860.0940.098
MODIS NDVI vs Actual OLI NDVIFused NDVI vs Actual OLI NDVI
Mean difference0.0240.032
Correlation coefficient (R)
p-value
0.8260.831
0.00010.0005

Share and Cite

MDPI and ACS Style

Zhang, M.; Lin, H.; Sun, H.; Cai, Y. RETRACTED: Estimation of Vegetation Productivity Using a Landsat 8 Time Series in a Heavily Urbanized Area, Central China. Remote Sens. 2019, 11, 133. https://doi.org/10.3390/rs11020133

AMA Style

Zhang M, Lin H, Sun H, Cai Y. RETRACTED: Estimation of Vegetation Productivity Using a Landsat 8 Time Series in a Heavily Urbanized Area, Central China. Remote Sensing. 2019; 11(2):133. https://doi.org/10.3390/rs11020133

Chicago/Turabian Style

Zhang, Meng, Hui Lin, Hua Sun, and Yaotong Cai. 2019. "RETRACTED: Estimation of Vegetation Productivity Using a Landsat 8 Time Series in a Heavily Urbanized Area, Central China" Remote Sensing 11, no. 2: 133. https://doi.org/10.3390/rs11020133

APA Style

Zhang, M., Lin, H., Sun, H., & Cai, Y. (2019). RETRACTED: Estimation of Vegetation Productivity Using a Landsat 8 Time Series in a Heavily Urbanized Area, Central China. Remote Sensing, 11(2), 133. https://doi.org/10.3390/rs11020133

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop