Next Article in Journal
Water Vapour Assessment Using GNSS and Radiosondes over Polar Regions and Estimation of Climatological Trends from Long-Term Time Series Analysis
Next Article in Special Issue
Monitoring Oasis Cotton Fields Expansion in Arid Zones Using the Google Earth Engine: A Case Study in the Ogan-Kucha River Oasis, Xinjiang, China
Previous Article in Journal
Cross-Comparison of Global Surface Albedo Operational Products-MODIS, GLASS, and CGLS
Previous Article in Special Issue
Alfalfa Yield Prediction Using UAV-Based Hyperspectral Imagery and Ensemble Learning
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Rapid Model (COV_PSDI) for Winter Wheat Mapping in Fallow Rotation Area Using MODIS NDVI Time-Series Satellite Observations: The Case of the Heilonggang Region

1
School of Earth Sciences and Resources, China University of Geosciences, Beijing 100083, China
2
Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100094, China
3
Guangdong Provincial Key Laboratory of Soil and Groundwater Pollution Control, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China
4
Institute at Brown for Environment and Society, Brown University, Providence, RI 02912, USA
*
Author to whom correspondence should be addressed.
Remote Sens. 2021, 13(23), 4870; https://doi.org/10.3390/rs13234870
Submission received: 14 October 2021 / Revised: 18 November 2021 / Accepted: 25 November 2021 / Published: 30 November 2021
(This article belongs to the Special Issue Smart Farming and Land Management Enabled by Remotely Sensed Big Data)

Abstract

:
Rapid and accurate monitoring of spatial distribution patterns of winter wheat over a long period is of great significance for crop yield prediction and farmland water consumption estimation. However, weather conditions and relatively long revisit cycles often result in an insufficient number of continuous medium-high resolution images over large areas for many years. In addition, the cropland pattern changes frequently in the fallow rotation area. A novel rapid mapping model for winter wheat based on the normalized difference vegetation index (NDVI) time-series coefficient of variation (NDVI_COVfp) and peak-slope difference index (PSDI) is proposed in this study. NDVI_COVfp uses the time-series index volatility to distinguish cultivated land from background land-cover types. PSDI combines the key growth stages of winter wheat phenology and special bimodal characteristics, substantially reducing the impact of abandoned land and other crops. Taking the Heilonggang as an example, this study carried out a rapid mapping of winter wheat for four consecutive years (2014–2017), and compared the proposed COV_PSDI with two state-of-the-art methods and traditional methods (the Spectral Angle Mapping (SAM) and the Iterative Self-Organizing Data Analysis Technique Algorithm (ISODATA)). The verification results revealed that the COV_PSDI model improved the overall accuracy (94.10%) by 4% compared with the two state-of-art methods (90.80%, 89.00%) and two traditional methods (90.70%, 87.70%). User accuracy was the highest, which was 93.74%. Compared with the other four methods, the percentage error (PE) of COV_PSDI for four years was the lowest in the same year, with the minimum variation range of PE being 1.6–3.6%. The other methods resulted in serious overestimation. This demonstrated the effectiveness and stability of the method proposed in the rapid and accurate extraction of winter wheat in a large area of fallow crop rotation region. Our study provides insight for remote sensing monitoring of spatiotemporal patterns of winter wheat and evaluation of “fallow rotation” policy implementation.

1. Introduction

Winter wheat is one of the most important food crops in the world and is directly related to global and regional food security issues [1,2]. In recent years, winter wheat cultivation is changing in some countries due to climate change, cropland conversion and crop planting structure adjustment [3]. The seeded acreage of winter wheat is highly influenced by fluctuations in food prices and agricultural markets [4]. Monitoring crop spatial patterns is essential for crop production management, investigating the contribution of agroecosystems to the terrestrial carbon cycle, and assessing the impact of global changes on regional agricultural production [5]. Winter wheat has a long growth cycle, with the first year planted in autumn and the second year harvested in summer. In the arid and semi-arid regions of the Northern Hemisphere, winter wheat is mainly grown under irrigation. The monoculture cropping pattern has led to severe over-exploitation of the groundwater table, which is particularly severe in northern China [6,7]. In response to this issue, China has adopted a policy of “fallow crop rotation” to reduce groundwater extraction and agricultural non-point source pollution [8,9]. Therefore, it is important to quickly and accurately obtain the spatiotemporal patterns of winter wheat during a long time for agricultural water consumption estimation, ecological function assessment and agricultural policy formulation [10,11].
Field surveys, statistical reporting and remote sensing monitoring are the main methods used to capture the accurate crop area and distribution patterns. Among them, remote sensing technology is currently one of the most effective and widely used monitoring methods due to its advantages of short cycle time, low cost and wide coverage [12,13,14,15]. Gilbertson et al. [16] applied the pan-sharpening multi-temporal Landsat 8 imagery to extract seven major crop types in the Cape Winelands region of South Africa by decision trees, k-nearest neighbor, support vector machine (SVM), and random forests (RF) methods. The results showed that the classification accuracy of the pan-sharpening was improved by 15%, and that SVM had consistently produced superior results. Jin et al. [17] used the normalized difference vegetation index (NDVI) time series data, which was derived from the Chinese HJ-1A/B satellite covering the whole growth period of winter wheat, to establish a classification model combined with time-integrated NDVI to map the irrigated and rainfed wheat areas. The application methods of remote sensing technology to obtain spatial and temporal distribution patterns of crops are mainly divided into three categories: thematic information extraction methods, vegetation index extraction methods and the comprehensive extraction methods of vegetation index and phenological information [18,19].
Thematic information extraction and vegetation index extraction methods, mainly based on single images of key crop growth periods, using support vector machine (SVM) [20], random forest (RF) [21] and decision tree [22] methods are used to carry out crop mapping at a field scale. Ideally, medium to high spatial resolution imagery can identify crop growth information with relative accuracy [23]. However, due to factors, such as weather conditions or satellite revisit cycles, it is often not possible to obtain continuous medium-high resolution remote sensing images over large areas for many years [24,25,26,27,28]. In addition, crop growth often changes significantly in the short term, and the selection of sample points and parameters on a single phase may lead to large uncertainties in terms of mapping accuracy [2]. The Moderate Resolution Imaging Spectroradiometer (MODIS) data can be effectively applied to the mapping of long-term and large-scale crops by employing its advantages of high temporal and spectral resolution and intermediate spatial resolution [5,29]. The high temporal resolution of MODIS can obtain more detailed and robust crop phenology information and capture the characteristics of crop growth change in the short term [30,31], which meets the data requirements of our model.
A comprehensive extraction method of vegetation index and phenological information based on phenological information can compensate for the shortcomings of single image data by effectively combining the crop spectral and seasonal rhythm information provided by multi-temporal images. For example, Dong et al. [32] used a time-weighted dynamic time warping method to map winter wheat by comparing the similarity of seasonal changes in the NDVI between different crops and standard winter wheat. Qu et al. [33] applied the NDVI time series data from the Landsat-8 operational land imager (OLI) to build the winter wheat index (WWI) based on the characteristics of changes around the overwintering and tassel stages, which improved the mapping accuracy. Sarvia et al. [34] applied the Minimum Distance and Random Forest classifiers to achieve a simple and effective mapping of wheat, with the NDVI multi-time series from Sentinel-2 satellite data. Time-series image extraction methods based on phenological information can quantify critical phenological features of the winter wheat growth cycle and eliminate the problem of crop phase mixing to a certain extent [17,29], but still with some limitations to some degree. On the one hand, studies which focus on the rapid and accurate monitoring of spatiotemporal patterns of winter wheat for large regional fallow and abandoned areas are insufficiently explored. On the other hand, existing models are generally poorly transplantable, have less consideration of the variability of crop phenology within large regions and the fluctuating characteristics of the temporal indices [35]. In addition, these studies have mostly assumed that the distribution of cropland is in a stable state for a long time, relying on the static land cover data to mask non-cropland information and reflect crop areas [36]. However, the distribution of cropland in the fallow rotation area changes frequently, resulting in the superposition of errors.
In light of these considerations, this study proposes a novel model that combines the coefficient of variation (COV) of NDVI with the bimodal characteristics of winter wheat for rapid extracting the spatiotemporal pattern of large-area winter wheat for fallow rotation areas. The specific contents include: (1) distinguishing cultivated land from background land-cover types by the sensitivity of COV to the fluctuation of time-series profile; (2) establishing the winter wheat rapid extraction index combining the crucial phenological information; (3) given the mixed effects of abandoned land and other crops, simulating the winter wheat in mixed pixels to determine the reasonable parameter threshold. A typical fallow rotation experimental area, Heilonggang, China, was selected as the study area to carry out large-scale rapid and accurate mapping of winter wheat for four consecutive years (2014–2017), and to evaluate the effectiveness of the model.

2. Materials and Methods

2.1. Study Area

The Heilonggang region is the main production area of winter wheat and one of the most important grain production bases in China [37]. It is located in the east-central part of the North China Plain, in the south-central region of Hebei Province (Figure 1), with an area of about 84,800 km2 and a population of about 51.45 million. The region is located in the warm temperate zone, with a continental monsoon climate and abundant light and heat resources. The traditional agricultural cropping pattern is mainly winter wheat-summer corn in a biannual system, as well as a few cotton and peanut, soybean-vegetable cropping patterns [3,38]. In addition, there are also some large areas of fruit trees, such as apricot trees, pear trees and so on. According to the field surveys and relevant phenological information (http://data.cma.cn/, accessed on 5 March 2020), the phenological periods of the main crops in the study area are summarized in Figure 2. Among them, winter wheat is the main overwintering crop in the study area [4,31,39].
Due to the monoculture cropping pattern, low-lying terrain and poor water drainage, a series of water and fertilizer resource problems have emerged in this area, such as shortage of water resources, deterioration of soil quality and deterioration of farmland ecological environment [40,41,42]. In recent years, with the implementation of the national “fallow rotation” policy, Heilonggang, a typical area in the groundwater leakage zone, has been used as a pilot area for the application of the fallow rotation system to carry out seasonal fallowing. The seasonal fallow is a “one season fallow, one season rainfed” system, whereby winter wheat, which needs to be pumped for irrigation, is left to fallow in order to reduce groundwater consumption and achieve sustainable agricultural development in the Heilonggang Plain [43].

2.2. Workflow

The detailed workflow for this study is schematized in Figure 3. First, we reconstructed the MODIS NDVI time-series dataset by the Savitzky-Golay (S-G) filter. Next, a rapid winter wheat extraction model (COV_PSDI) was constructed based on the phenological information analysis and NDVI time-series curve characteristics. The sensitivity of COV to NDVI curve fluctuation was used to separate cropland and non-cropland. The difference method was adopted to extract the bimodal characteristics of winter wheat. Following this, combined with the constructed slope difference index (SDI), the winter wheat in the fallow rotation area was accurately extracted. Third, we selected the sample points through the preprocessed high-resolution images. The reasonable thresholds for the parameters were determined by simulating the proportion of winter wheat in mixed pixels. Finally, the accuracy of the experimental results were evaluated and compared by calculating the overall accuracy (OA), producer accuracy (PA), user accuracy (UA), and percentage error (PE) of different methods, to verify the feasibility of the method.

2.3. Data and Data Pre-Processing

2.3.1. Data

NDVI data from the MODIS MOD13Q1 product were used for the crop growth curves analysis and winter wheat mapping in this study. MOD13Q1 product were synthesized over 16 days based on the maximum value composite (MVC) method, and were obtained from the U.S. National Aeronautics and Space Administration (NASA) website (https://ladsweb.modaps.eosdis.nasa.gov/search/, accessed on 25 April 2020) (Table 1). The dataset was corrected geometrically and atmospherically. This NDVI data captures vegetation phenology information, which has been widely used in crop mapping [44,45,46].
Landsat-8 OLI multi-spectral image, Google Earth images and statistics data were used to sample points selection and model validation. Landsat-8 OLI images with 30 m spatial resolution were downloaded from the U.S. Geological Survey (USGS) website (http://earthexplorer.usgs.gov/, accessed on 10 April 2020). The acquisition dates of the four scenes of Landsat-8 OLI images used in this study covered three important growth nodes of winter wheat, respectively. Specifically, December 2016 or January 2017 (overwintering stages), April or May 2017 (heading stages) and June 2017 (harvesting stages). Statistical data of winter wheat area were derived from the statistical yearbooks of Hebei Province and cities from 2014 to 2017 (http://tjj.hebei.gov.cn/, accessed on 26 July 2020).
Sample regions were selected from OLI images. The sample areas, i.e., the yellow areas in Figure 4, were delineated according to the image quality (the three images of three key phenological periods have no cloud cover). The sample pixels accounted for about 25% of the OLI images (Figure 4). Considering the resolution difference between MODIS NDVI products and validation dataset, sample points were selected at the center of at least 16 × 16 pure pixels in Landsat-8 OLI images to ensure that the sample points were pure pixels in MODIS NDVI images. Based on a field survey, Google Earth history images, and OLI images of three key phenological periods (Figure 4i–iii), the sample points were selected in each yellow region. Two groups of independent reference samples were selected. One was used to analyze the NDVI time-series curve of different land cover types and threshold simulation analysis to develop a rapid extraction model of winter wheat. The other was used to validate the winter wheat mapping results. The reference sample set for verification was randomly selected in the delineated yellow areas. The reference samples of each pixel were visually labeled as winter wheat or non-winter wheat. Finally, 492 winter wheat points and 508 non-winter wheat points were used to obtain mapping accuracy.

2.3.2. Data Preprocessing

MODIS Reprojection Tool (MRT) was used to batch process the NDVI time series data. Image mosaic, reprojection (the coordinate system was WGS84, and the projection system was UTM) and clipping were performed on the 4-track data covering the study area. To improve the spatial coregistration between the two sensors, geographic registration was performed in ENVI 5.3. In addition, the S-G filter was used to smooth the NDVI time series curve [48,49]. This filtering method can effectively remove the influence of cloud, atmosphere and data quality in MODIS images. Take winter wheat-summer corn-fallow-summer corn sample points as an example analysis (Figure 5). The reconstructed NDVI time series curve has eliminated jagged irregular fluctuations, and clearly described the changing trend characteristics and local mutation information of the NDVI time series. Moreover, it conformed to the growth characteristics of crops and met the requirements of identifying different crop planting patterns. After all of the preprocessing steps were completed, the processed MODIS NDVI time series dataset was used for subsequent processing and analysis.

2.4. NDVI Time Series of Different Land Cover Types

NDVI is widely regarded as a vegetation index that is able to effectively retrieve vegetation information [34,50], especially in phenology, tree vigor decline evaluation, and crop-yield forecasting [51,52,53]. NDVI was selected as the vegetation phenology predictor in this study. The mean NDVI time series profiles of different land cover types were produced using the selected reference samples. Figure 6 shows that the cropland NDVI time series curve varies obviously with the cultivation and harvest of crops, especially the two-season planting pattern (i.e., winter wheat-summer corn, minor crops-vegetables), which are also found in existing studies [33,54].
The temporal NDVI of winter wheat-summer corn contained two large peaks and one small peak. The first big peak (i.e., PI1 in Figure 6) appeared from the heading stage to the filling stage, while other crops and natural vegetation were in the growth stage. The other vegetation NDVI showed an upward trend, and was lower than winter wheat. Subsequently, the winter wheat entered the mature stage, and the NDVI curve dropped sharply. In addition to the fallow farmland, the NDVI of other crops and natural vegetation were still in the rising stage to varying degrees. The fallow farmland NDVI increased due to the growth of weeds [32]. At the early stage of summer corn sowing, weed eradication caused a significant decline in the NDVI curve, but the slope of the winter wheat curve was higher than that of weeds. The second large peak occurred at phase 15 (July), when the NDVI of summer corn went to the highest value. During this period, the cotton and forest NDVI values were also high, but there were no peaks. At the phase 22 (overwintering period), the winter wheat wavelet peak PI2 appeared. At this time, most other crops were harvested, the leaves of the forest were withered, and the NDVI value of natural vegetation was low [33].
The NDVI time series curve of small crops (soybean and peanut) and vegetables (cabbage and broccoli) contained two peaks, which was also based on a two-season planting mode in one year. However, due to the different phenological periods, the NDVI peak was significantly different from that of the winter wheat-summer corn [5]. The NDVI curves of single-season cotton and fallow farmland-summer corn had obvious phenological characteristics of the crop planting pattern with one seasonal crop a year, while the NDVI curves only contained one peak. There were very few winter wheat-oil sunflower/soybean planting areas in the study area, as these were combined with winter wheat-summer corn in this study.
The NDVI time series curve of forests (including orchards) began to increase rapidly after phase 6, and the NDVI value was approximately 0.8 at phase 10 to 16 and remained stable. After phase 17, the leaves withered, and the NDVI began to decline slowly. The NDVI curve of construction land (including bare land) was low, with no more than 0.3, and there was no obvious peak valley phenomenon.

2.5. Winter Wheat Extraction Model (COV_PSDI)

2.5.1. Variation Coefficient of Fluctuation Period NDVI_COVfp

NDVI_COV represents the variation coefficient of the NDVI time series, that is, the ratio of the standard deviation to mean of the NDVI time-series of each pixel [55,56]. NDVI_COV is used to eliminate non-cultivated land information and reduce the interference of background land cover information, such as woodland (orchard), grassland and construction land. As shown in Figure 6, the NDVI curve of construction land is at a low value throughout the year. The growth of forest and grassland is stable from May to September (i.e., phase 9 to 16), and the NDVI time series curve fluctuation is small. However, during this period, due to crop harvesting and sowing, the NDVI time series curve of cultivated land fluctuates greatly. This study defines this period as a fluctuation period. Therefore, the NDVI_COVfp is constructed by using the sensitivity of COV to NDVI curve fluctuation.
σ = 1 n 1 [ k = 1 n ( N D V I k ) 2 1 n ( k = 1 n N D V I k ) 2 ] ,
μ = k = 1 n ( N D V I k ) / n ,
where σ and μ are the standard deviation and mean value of NDVI of each pixel in phase 9 to 16 (i.e., DOY129 to DOY241), respectively; NDVIk are NDVI values of phase 9 to 16 (i.e., DOY129 to DOY241); and n is the number of phases. The COV calculation formula of the pixel for row i and column j is given by:
C O V i j = σ i j / μ i j ,
The smaller the COV value of the NDVI time series, the smoother the fluctuation of the NDVI time series curve.

2.5.2. Peak and Slope Difference Index PSDI

As described in Section 2.4, as the typical overwintering crop in the study area, the NDVI time series curve of winter wheat shows special bimodal characteristics over the whole growth period. The latitude difference in the study area is small, which can be considered as being because the growth period of winter wheat is basically the same. NDVI peaks at tillering and heading-filling stages of winter wheat are extracted by a difference method, which is a method used to obtain the maximal value of discrete points [57]:
U = N D V I t N D V I t + 1 ,
U > 0 , V = 1 ; U < 0 , V = 1 ,
P I = V t V t + 1 ,
where U is the NDVI difference between the NDVI value of phase t and phase t + 1; V = 1 means the NDVI of phase t is greater than that of phase t + 1, and V = −1 means the NDVI of phase t is less than that of phase t + 1; PI = 2 is the minimal of the time series curve (trough), and PI = −2 is the maximal of the time series curve (peak).
It should be mentioned that the quadratic difference is very sensitive to extreme values. It can extract each small extreme value [57]. The slope difference index (SDI in Figure 6) is constructed to avoid the occurrence of error extremum, which has great potential for identifying winter wheat. Considering the inter-annual differences in crop growth and climatic factors that can lead to delayed and advanced phenology, NDVIi and NDVIj are dynamically selected for the maximum value in the three adjacent phases. The SDI is defined as:
S D I = N D V I j N D V I i N D V I i ,
{ N D V I i = m a x { N D V I t | i 1 t i + 1 } N D V I j = m i n { N D V I t | j 1 t j + 1 } ,
where NDVIi is the inflection point value of NDVI curve at winter wheat heading stage; NDVIj is the NDVI value after the harvest.

2.6. Threshold Analysis

The reasonable setting of NDVI_COVfp and SDI threshold is crucial for winter wheat mapping [58]. Compared with woodland, construction land and bare land, the NDVI curve of winter wheat–summer corn fluctuates greatly during the fluctuation period, and the NDVI_COVfp value (0.34) is high (Figure 7). However, the NDVI_COVfp value of small crops–vegetables is also high, and winter wheat could not be effectively extracted only through the NDVI_COVfp. The SDI value of winter wheat is negative and lower than −0.50. Combining the NDVI_COVfp and SDI with the bimodal characteristic of winter wheat, the extraction of pure winter wheat pixels is relatively simple. However, selecting a reasonable threshold to extract winter wheat from mixed pixels is a difficult issue. In order to simplify the mixed pixel issue, this study defines that, if most of the area (50% or more) in a pixel (250 m × 250 m) is covered by winter wheat, it is defined as a pure winter wheat pixel [59]. Based on the NDVI values of pure pixels of different land cover types, the threshold is determined by analyzing the index values and the percentage of winter wheat in the mixed pixels with 10% winter wheat mixing ratio as the increment. Whereas, the NDVI values of mixed pixels with various area percentages were detected by simulation analysis based on the NDVI values of pure pixels of different land cover types. For instance, assuming that the NDVI values of winter wheat pure pixels and forest pure pixels in phase t are NDVIWheat,t and NDVIForest,t, respectively. The NDVI values of mixed pixels with 90% of winter wheat and 10% of forest in phase t are NDVImixed = 0.9 × NDVIWheat,t + 0.1 × NDVIForest,t. Afterward, the parameter values were calculated based on the simulated NDVI values of various mixed pixels in different periods.

2.7. Accuracy Assessment and Model Comparison

The proposed COV_PSDI is evaluated in two aspects, i.e., qualitative evaluation and quantitative evaluation. The confusion matrix and percentage error are used as verification standards [60,61,62]. As a standard format for the accuracy evaluation of remote sensing classification models, the confusion matrix includes OA, PA, and UA [63,64,65].
For quantitative evaluation, images extraction results and statistical yearbook data are compared using the percentage error method to assess the proposed model effectiveness.
P E = | E A r e S A r e | S A r e × 100 % ,
where EAre is the winter wheat area based on different winter wheat mapping methods in this study; and SAre is the winter wheat area recorded in the statistical yearbook.
To validate the effectiveness of the proposed COV_PSDI, the proposed method is compared with two existing winter wheat mapping methods using different key characteristics of the NDVI time series profiles (Table 2). In this study, the NDVI time-series data were used instead of enhanced vegetation index (EVI) data in the CBAH method. Thus, two indices, the NDVI change amplitude during seedling stage to heading stage and NDVI change amplitude during heading stage to harvesting stage were used. In addition, two conventional image classification methods (i.e., the Spectral Angle Mapping (SAM) [66,67] and Iterative Self-Organizing Data Analysis Technique Algorithm (ISODATA) [68,69]) are applied to winter wheat mapping for comparison. SAM achieves rapid mapping by calculating the similarity between the endmember spectrum and the pixel spectrum. The smaller angle between the two spectra, the higher the similarity [70,71]. Due to the SAM algorithm utilizing only the vector direction rather than the vector length, it performs well in regions of uniform surface coverage type. In this paper, the NDVI time-series profile was used as a spectrum to permit the rapid mapping of winter wheat. ISODATA is one of the commonly used unsupervised classification methods. It allows specifying different numbers of clusters (i.e., ranging from the minimum to maximum number of clusters) [72], which is more adaptable and flexible than K-means.

3. Results

3.1. Analysis Results of Model Parameters

The relationships between NDVI_COVfp, SDI, WWI, CBAH, and the proportion of winter wheat in mixed pixels are displayed in Figure 8. When winter wheat was mixed with non-farmland, such as forest, grassland, bare land and construction land, the NDVI_COVfp values decreased gradually with the decrease of the winter wheat proportion in mixed pixels. On the contrary, the SDI values increased with the decrease of the winter wheat proportion in mixed pixels. When the proportion of winter wheat was over 50%, the values of NDVI_COVfp and SDI varied slowly. After the proportion of winter wheat was less than 50%, the NDVI_COVfp and SDI values decreased and increased rapidly, respectively (Figure 8c,d). Combined with the actual investigation in the study area (the mixed pixel is mainly composed of winter wheat and other crops, winter wheat and forest), the model NDVI_COVfp threshold value should be 0.21, and the SDI index value should be −0.16.
The same threshold analysis method as NDVI_COVfp and SDI was used to analyze the thresholds of the CBAH and WWI models. The CBAH values were the sum of EVE and EVL in Figure 8c. This dropped linearly with the decrease of the winter wheat proportion in mixed pixels, with the maximum variation. According to the definition of mixed pixels in Section 3.3 and the actual investigation in the study area (the mixed pixel is mainly composed of winter wheat and other crops, winter wheat and forest), the model CBAH threshold value should be 0.34 in the simulation analysis results. The WWI values showed small differences in different types of winter wheat mixed pixels (Figure 8d). However, in the mixed pixels of winter wheat and cotton, the WWI value decreased rapidly after the proportion of winter wheat was less than 50%. The WWI values of other winter wheat mixed pixels were stable. Therefore, WWI = 0.005 was selected in this paper.

3.2. Winter Wheat Mapping Results

Based on the COV_PSDI method and threshold analysis, the spatial mapping of winter wheat in Heilonggang and its surrounding areas from 2014 to 2017 was completed (Figure 9). All pixels were classified into two categories, i.e., winter wheat and non-winter wheat. Figure 9 illustrated that the spatial distribution of winter wheat was quite different. Winter wheat is mainly concentrated in the central and western plains, Baoding, Shijiazhuang, Handan, central Xingtai, southern Cangzhou and Hengshui. The reasons for the lower distribution of winter wheat in the west and northeast might be that the west of Heilonggang is mostly mountains, and land salinization is a serious issue in northern coastal areas of Cangzhou. Therefore, there are few suitable areas for winter wheat planting. This distribution pattern is consistent with the existing literature [36,73].
In 2014 and 2015, the total area of winter wheat in Heilonggang remained basically stable, being 22,076.2 km2 and 22,072.8 km2 in terms of proposed COV_PSDI models, respectively. With the implementation of the “fallow rotation” policy, the planting information of winter wheat in Hengshui and Cangzhou changed greatly, and the area decreased significantly. The total area of winter wheat in the Heilonggang area decreased to 20,416.8 km2 in 2016. However, owing to the long-term traditional planting pattern in the study area, as well as the insufficient subsidy for fallow or the rise in food prices [74], some farmers did not fallow and the area of winter wheat increased in 2017. The pilot of the fallow project in Heilonggang should be continued and monitored regularly.

3.3. Model Comparison and Validation

The mapping accuracy of the five models in 2017 was verified with the sample points from Landsat 8 and Google Earth images (Table 3). For the 492 winter wheat sample points, 464 points were correctly classified (94.32% agreement) in the COV_PSDI result. The mapping accuracy of our model was improved by 4% on OA, compared with the two state-of-art methods (CBAH: 90.80%, WWI: 89.00%) and two traditional methods (SAM: 90.70%, ISOData: 87.70%). The UA of the COV_PSDI were 93.74%, which were significantly higher than those of the other four models.
Although the PA of the COV_PSDI model (94.32%) was lower than that of the WWI model (99.39%), it was higher than that of the SAM (91.67%), ISODATA (92.89%), and CBAH (93.70%) models. The UA of the WWI model was 82.05%, indicating that many non-winter wheat pixels were misclassified as winter wheat (i.e., 107 points were misclassified in the 508 non-winter wheat sample points), which greatly reduced the OA of the WWI model. It fully shows that our method had certain advantages in winter wheat mapping in the fallow area.
We also compared the mapping results of the different models with the study area’s statistical data concerning winter wheat. Table 4 showed that a slight overestimation from COV_PSDI was revealed when compared with the areas of winter wheat from the statistical data in 2014, 2015, and 2017. Underestimates were found in the COV_PSDI results for 2016. However, SVM, ISODATA and WWI seriously overestimated the results, especially ISODATA, with PE reached to 37.3% in 2015. On the contrary, there were greatly underestimates in the CBAH results, with PE as high as 25.2% in 2015. The PE of COV_PSDI winter wheat mapping results were 2.1%, 3.6%, 1.6% and 3.4% from 2014 to 2017, respectively. Compared with other models, the PE of the model proposed in this study was the lowest in the same year, and the variation range was the smallest (1.6–3.6%). The changes of PE in the other models were SAM: 3.9–10.9%, ISODATA: 4.3–37.3%, CBAH: 3.0–25.2%, and WWI: 7.0–20.3%.
Overall, the winter wheat mapping results of the COV_PSDI model were more consistent with the actual situation. The COV_PSDI model was more stable in the mapping research of winter wheat for large-scale fallow rotation areas.

4. Discussion

The time series classification model considers the key phenology characteristics of vegetation [75,76]. The proposed COV_PSDI model achieved rapid winter wheat mapping in fallow rotation areas with frequent changes of cropland distribution and improved accuracy. As shown in our study, COV_PSDI outperformed the existing time-series indices (CBAH and WWI) [33,36] for winter wheat and mapping methods (SAM and ISODATA) [77,78]. In CBAH, the NDVI values of the three growth stages (seeding, heading and harvesting) were utilized for winter wheat mapping. The agricultural land in the North China Plain is mainly irrigated [17]. Therefore, watering in April of the orchard would make the NDVI time-series curve drop suddenly during phase 8 and 9 (i.e., the heading stage of winter wheat), resulting in a misclassification of the CBAH index (Figure 10e). The ISODATA method also has this misclassification phenomenon (Figure 10d), owing to it only considering the shape of the time-series curve [79]. The SAM and ISODATA methods (Table 4) have low adaptability [80]. Large area or long-term mapping needs to select time-series curve sample points or parameters alone. In the study of Nguyena et al. [81], the accuracy of identifying cropland was improved by submitting phenological indices to the random forest classifier, but the accuracy of identifying wheat was only 93%. In addition, because of frequent changes of cropland in the fallow rotation area, cumulative errors due to relying on the static land cover data to mask non-cropland information would increase the uncertainties associated with the mapping outcomes. Wang et al. [82] extracted the multiple cropping index from 2006 to 2011 based on the 2010 cropland boundary data obtained by the mask, resulting in the failure to detect the multiple types of cropping information in the fallow areas, limiting the extraction accuracy. Our model, COV_PSDI, avoided the above misclassification, which not only improved the mapping accuracy of winter wheat, but also ensured the stability of multi-year extraction results.
COV is sensitive to curve variation, and the selection of the fluctuation period is very important in order to eliminate the interference of non-farmland background information. It can be seen from this study that the annual NDVI time-series curve of forest and grassland also experienced great fluctuation (Figure 6). We calculated the annual variation coefficient of NDVI time-series profiles (NDVI_COVannual) of different land cover types (Figure 11). It is obvious that the variance of NDVI_COVannual between cultivated land and non-cultivated land is small, which reduces the sensitivity of NDVI_COV to cropland identification, whereas the NDVI_COVfp calculated with May to September as the fluctuation period can effectively distinguish cultivated land and non-cultivated land information (Figure 7).
Compared with EVI, NDVI was easily saturated at high vegetation coverage. However, NDVI can better characterize the surface vegetation coverage information in low vegetation coverage [83]. The winter wheat mapping model proposed in this paper mainly used several important phenology nodes of winter wheat, i.e., overwintering period (medium coverage), heading period (high coverage) and harvest period (low coverage). NDVI changed obviously before and after the heading period of winter wheat (Figure 6). Therefore, the model was rarely affected by NDVI saturation. Moreover, NDVI is the most sensitive to medium vegetation coverage, and COV is directly affected by the change of vegetation coverage [55]. Some mixed pixels with large proportion sparse woodland (orchard) or fallow land showed lower NDVI values during winter wheat harvest. Ignoring the COV constraint will lead to misclassification of winter wheat pixels, because the SDI value is greater than the threshold. Although the classification accuracy of COV_PSDI is the highest, the results are different for other regions with different vegetation coverage characteristics. This is because the sensitivity of NDVI to different vegetation cover types is different [84]. The fluctuation period of the winter wheat NDVI curve and the exact date of three key growth stages are different in different geographical locations, but the growth period of winter wheat itself and its corresponding NDVI characteristics are similar in different regions. Therefore, the COV_PSDI is also suitable for large-scale winter wheat mapping with significant regional differences.
Although the proposed COV_PSDI produced promising results, there are still several shortcomings. Firstly, due to the limitation of MODIS imaging quality, NDVI time series data are greatly affected by clouds. Whether or not cloud removal is carried out, there still be some thin clouds that are not completely eliminated, resulting in an NDVI value that is slightly lower than the real value. Second, the mixed pixel issue limits the mapping accuracy of winter wheat. Although the mixed pixel problem is considered in the threshold analysis, it still inevitably affects the extraction accuracy. Third, owing to technical limitations (mainly adopting the method of layer-by-layer reporting) and historical background, statistical yearbook data are underreported and misreported, which is another reason for inconsistent classification results.
In the future, integrating multi-source multispectral sensor data, such as Landsat 8/OLI and Sentinel-2/MSI, can be considered to produce fusion data with uniform high spatial resolution and high temporal resolution for many years [85,86,87,88], which will be expected to improve the accuracy of large area mapping. How to obtain large-scale and long-term continuous high resolution remote sensing images is also an urgent problem in this field. Our study only compares the exponential model method. Machine learning and deep learning methods can be further used to evaluate winter wheat mapping in fallow rotation areas. In addition, this study mapped winter wheat in typical fallow rotation areas based on the COV_PSDI method. However, under the influence of fallow rotation policy, how to carry out fallow and the effect of fallow have not been discussed, and these issues need to be further refined and clarified in future research.

5. Conclusions

In this study, we proposed a new approach, COV_PSDI, for the spatiotemporal pattern extraction of winter wheat for a large area of a fallow rotation area. The unique phenological period of winter wheat and the winter wheat-summer corn planting pattern are the theoretical basis for the remote sensing mapping of winter wheat based on an NDVI time-series dataset. The COV of the NDVI time series was applied to separate farmland from background land cover types, and the difference method was used to extract the bimodal characteristics of winter wheat. Given the sensitivity of the quadratic difference method to extreme value, SDI was constructed. One effective mapping study of winter wheat in the large-scale fallow rotation region in the Heilonggang area was carried out. The results demonstrated that the variation coefficient of the fluctuation period (NDVI_COVfp) showed better separability between cultivated land and non-cultivated land. The overall accuracy and user accuracy of the COV_PSDI were high, at 94.10% and 93.74%, respectively. Using statistical data to verify the model, COV_PSDI was slightly overestimated, with a variation range of PE of 1.6–3.6%. Moreover, the COV_PSDI model was more stable and accurate, compared with four comparative methods. Therefore, the proposed COV_PSDI model provides a novel method for mapping winter wheat in large-scale fallow rotation areas, and may be applicable to other study areas.

Author Contributions

Conceptualization, X.Z. and K.L.; methodology, X.Z. and S.W.; software, X.Z.; validation, X.Z., S.W. and K.L.; formal analysis, X.Z.; investigation, X.Z.; resources, S.W., K.L. and X.L. (Xin Long); data curation, S.W. and K.L.; writing—original draft preparation, X.Z.; writing—review and editing, X.Z., K.L. and X.L. (Xueke Li); visualization, X.Z.; supervision, X.L. (Xueke Li); project administration, X.Z.; funding acquisition, K.L., S.W and X.L. (Xin Long). All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the National Natural Science Foundation of China, grant number 41907191 and 41671362, the Shenzhen Science and Technology Innovation Committee, grant number JCYJ20190809144601730, and Special Fund for Transformation of Scientific and Technological Achievements in Inner Mongolia Autonomous Region, grant number 2021CG0045.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

A list of all abbreviations in alphabetical order to define all abbreviations used in this study.
AbbreviationDescription
CBAHcombining variations before and after estimated heading dates
CLCDannual China Land Cover Dataset
COV coefficient of variation
DOYthe day of a year
EVEearly growth stage (seeding stage to heading stage)
EVIenhanced vegetation index
EVLlate growth stage (heading stage to harvesting stage)
ISODATAiterative self-organizing data analysis technique algorithm
MODISModerate Resolution Imaging Spectroradiometer
MRTMODIS reprojection tool
MVCmaximum value composite
NASANational Aeronautics and Space Administration
NDVInormalized difference vegetation index
NDVI_COVvariation coefficient of time-series NDVI
NDVI_COVannualannual variation coefficient of NDVI time-series
NDVI_COVfpvariation coefficient of time-series NDVI in fluctuation period
NIRnear-infrared
OAoverall accuracy
OLIoperational land imager
PAproducer accuracy
PEpercentage error
PIpeak index
PSDIpeak-slope difference index
RFrandom forests
SAMspectral angle mapping
SDIslope difference index
S-G Savitzky-Golay
SVMsupport vector machine
UAuser accuracy
USGSU.S. Geological Survey
WWIwinter wheat index

References

  1. Chen, Y.; Zhang, Z.; Tao, F. Improving regional winter wheat yield estimation through assimilation of phenology and leaf area index from remote sensing data. Eur. J. Agron. 2018, 101, 163–173. [Google Scholar] [CrossRef]
  2. Zhang, D.; Fang, S.; She, B.; Zhang, H.; Jin, N.; Xia, H.; Yang, Y.; Ding, Y. Winter Wheat Mapping Based on Sentinel-2 Data in Heterogeneous Planting Conditions. Remote Sens. 2019, 11, 2647. [Google Scholar] [CrossRef] [Green Version]
  3. Zhou, T.; Pan, J.; Zhang, P.; Wei, S.; Han, T. Mapping Winter Wheat with Multi-Temporal SAR and Optical Images in an Urban Agricultural Region. Sensors 2017, 17, 1210. [Google Scholar] [CrossRef] [PubMed]
  4. Zhang, X.; Qiu, F.; Qin, F. Identification and mapping of winter wheat by integrating temporal change information and Kullback–Leibler divergence. Int. J. Appl. Earth. Obs. 2019, 76, 26–39. [Google Scholar] [CrossRef]
  5. Tao, J.; Wu, W.; Zhou, Y.; Wang, Y.; Jiang, Y. Mapping winter wheat using phenological feature of peak before winter on the North China Plain based on time-series MODIS data. J. Integr. Agric. 2017, 16, 348–359. [Google Scholar] [CrossRef] [Green Version]
  6. Wang, J.; Jiang, Y.; Wang, H.; Huang, Q.; Deng, H. Groundwater irrigation and management in northern China: Status, trends, and challenges. Int. J. Water Resour. Dev. 2020, 36, 670–696. [Google Scholar] [CrossRef]
  7. Wu, X.; Qi, Y.; Shen, Y.; Yang, W.; Zhang, Y.; Kondoh, A. Change of winter wheat planting area and its impacts on groundwater depletion in the North China Plain. J. Geogr. Sci. 2019, 29, 891–908. [Google Scholar] [CrossRef] [Green Version]
  8. Li, Q.; Ma, L.Y.; Wang, X.F.; Wei, M.; Shi, Q.H.; Yang, F.J. Effects of short-term summer fallow and rotations on soil quality under plastic greenhouse cultivation. J. Food Agric. Environ. 2012, 10, 1106–1110. [Google Scholar]
  9. Yang, X.; Steenhuis, T.; Davis, K.; Van der Werf, W.; Ritsema, C.; Pacenka, S.; Zhang, F.; Siddique, K.; Du, T. Diversified crop rotations enhance groundwater and economic sustainability of food production. Food. Energy Secur. 2021, 10, e311. [Google Scholar] [CrossRef]
  10. Skakun, S.; Franch, B.; Vermote, E.; Roger, J.-C.; Becker-Reshef, I.; Justice, C.; Kussul, N. Early season large-area winter crop mapping using MODIS NDVI data, growing degree days information and a Gaussian mixture model. Remote Sens. Environ. 2017, 195, 244–258. [Google Scholar] [CrossRef]
  11. Wu, D.; Yu, Q.; Lu, C.; Hengsdijk, H. Quantifying production potentials of winter wheat in the North China Plain. Eur. J. Agron. 2006, 24, 226–235. [Google Scholar] [CrossRef]
  12. Khan, A.; Hansen, M.C.; Potapov, P.V.; Adusei, B.; Pickens, A.; Krylov, A.; Stehman, S.V. Evaluating Landsat and RapidEye Data for Winter Wheat Mapping and Area Estimation in Punjab, Pakistan. Remote Sens. 2018, 10, 489. [Google Scholar] [CrossRef] [Green Version]
  13. Liu, J.; Zhu, W.; Atzberger, C.; Zhao, A.; Pan, Y.; Huang, X. A Phenology-Based Method to Map Cropping Patterns under a Wheat-Maize Rotation Using Remotely Sensed Time-Series Data. Remote Sens. 2018, 10, 1203. [Google Scholar] [CrossRef] [Green Version]
  14. Wu, M.; Huang, W.; Niu, Z.; Wang, Y.; Wang, C.; Li, W.; Hao, P.; Yu, B. Fine crop mapping by combining high spectral and high spatial resolution remote sensing data in complex heterogeneous areas. Comput. Electron. Agric. 2017, 139, 1–9. [Google Scholar] [CrossRef]
  15. Li, X.; Wu, T.; Liu, K.; Li, Y.; Zhang, L. Evaluation of the Chinese Fine Spatial Resolution Hyperspectral Satellite TianGong-1 in Urban Land-Cover Classification. Remote Sens. 2016, 8, 438. [Google Scholar] [CrossRef] [Green Version]
  16. Gilbertson, J.K.; Kemp, J.; van Niekerk, A. Effect of pan-sharpening multi-temporal Landsat 8 imagery for crop type differentiation using different classification techniques. Comput. Electron. Agric. 2017, 134, 151–159. [Google Scholar] [CrossRef] [Green Version]
  17. Jin, N.; Tao, B.; Ren, W.; Feng, M.; Sun, R.; He, L.; Zhuang, W.; Yu, Q. Mapping Irrigated and Rainfed Wheat Areas Using Multi-Temporal Satellite Data. Remote Sens. 2016, 8, 207. [Google Scholar] [CrossRef] [Green Version]
  18. Qin, Y.; Xiao, X.; Dong, J.; Zhou, Y.; Zhu, Z.; Zhang, G.; Du, G.; Jin, C.; Kou, W.; Wang, J.; et al. Mapping paddy rice planting area in cold temperate climate region through analysis of time series Landsat 8 (OLI), Landsat 7 (ETM+) and MODIS imagery. ISPRS J. Photogramm. 2015, 105, 220–233. [Google Scholar] [CrossRef] [Green Version]
  19. Fang, H.; Wei, Y.; Dai, Q. A Novel Remote Sensing Index for Extracting Impervious Surface Distribution from Landsat 8 OLI Imagery. Appl. Sci. 2019, 9, 2631. [Google Scholar] [CrossRef] [Green Version]
  20. Zheng, B.; Myint, S.W.; Thenkabail, P.S.; Aggarwal, R.M. A support vector machine to identify irrigated crop types using time-series Landsat NDVI data. Int. J. Appl. Earth Obs. 2015, 34, 103–112. [Google Scholar] [CrossRef]
  21. Zhong, L.; Gong, P.; Biging, G.S. Efficient corn and soybean mapping with temporal extendability: A multi-year experiment using Landsat imagery. Remote Sens. Environ. 2014, 140, 1–13. [Google Scholar] [CrossRef]
  22. Tian, H.; Huang, N.; Niu, Z.; Qin, Y.; Pei, J.; Wang, J. Mapping Winter Crops in China with Multi-Source Satellite Imagery and Phenology-Based Algorithm. Remote Sens. 2019, 11, 820. [Google Scholar] [CrossRef] [Green Version]
  23. Martínez-Casasnovas, J.A.; Martín-Montero, A.; Auxiliadora Casterad, M. Mapping multi-year cropping patterns in small irrigation districts from time-series analysis of Landsat TM images. Eur. J. Agron. 2005, 23, 159–169. [Google Scholar] [CrossRef]
  24. Castillejo-González, I.L.; López-Granados, F.; García-Ferrer, A.; Peña-Barragán, J.M.; Jurado-Expósito, M.; de la Orden, M.S.; González-Audicana, M. Object- and pixel-based analysis for mapping crops and their agro-environmental associated measures using QuickBird imagery. Comput. Electron. Agric. 2009, 68, 207–215. [Google Scholar] [CrossRef]
  25. Conese, C.; Maselli, F. Use of multitemporal information to improve classification performance of TM scenes in complex terrain. ISPRS J. Photogramm. 1991, 46, 187–197. [Google Scholar] [CrossRef]
  26. Whitcraft, A.K.; Vermote, E.F.; Becker-Reshef, I.; Justice, C.O. Cloud cover throughout the agricultural growing season: Impacts on passive optical earth observations. Remote Sens. Environ. 2015, 156, 438–447. [Google Scholar] [CrossRef]
  27. Liu, K.; Su, H.; Zhang, L.; Yang, H.; Zhang, R.; Li, X. Analysis of the Urban Heat Island Effect in Shijiazhuang, China Using Satellite and Airborne Data. Remote Sens. 2015, 7, 4804. [Google Scholar] [CrossRef] [Green Version]
  28. Li, X.; Liu, K.; Tian, J. Variability, predictability, and uncertainty in global aerosols inferred from gap-filled satellite observations and an econometric modeling approach. Remote Sens. Environ. 2021, 261, 112501. [Google Scholar] [CrossRef]
  29. Wardlow, B.D.; Egbert, S.L.; Kastens, J.H. Analysis of time-series MODIS 250 m vegetation index data for crop classification in the U.S. Central Great Plains. Remote Sens. Environ. 2007, 108, 290–310. [Google Scholar] [CrossRef] [Green Version]
  30. Becker-Reshef, I.; Vermote, E.; Lindeman, M.; Justice, C. A generalized regression-based model for forecasting winter wheat yields in Kansas and Ukraine using MODIS data. Remote Sens. Environ. 2010, 114, 1312–1323. [Google Scholar] [CrossRef]
  31. Pan, Y.; Li, L.; Zhang, J.; Liang, S.; Zhu, X.; Sulla-Menashe, D. Winter wheat area estimation from MODIS-EVI time series data using the Crop Proportion Phenology Index. Remote Sens. Environ. 2012, 119, 232–242. [Google Scholar] [CrossRef]
  32. Dong, J.; Fu, Y.; Wang, J.; Tian, H.; Fu, S.; Niu, Z.; Han, W.; Zheng, Y.; Huang, J.; Yuan, W. Early-season mapping of winter wheat in China based on Landsat and Sentinel images. Earth Syst. Sci. Data. 2020, 12, 3081–3095. [Google Scholar] [CrossRef]
  33. Qu, C.; Li, P.; Zhang, C. A spectral index for winter wheat mapping using multi-temporal Landsat NDVI data of key growth stages. ISPRS J. Photogramm. 2021, 175, 431–447. [Google Scholar] [CrossRef]
  34. Sarvia, F.; Xausa, E.; De Petris, S.; Cantamessa, G.; Borgogno-Mondino, E. A Possible Role of Copernicus Sentinel-2 Data to Support Common Agricultural Policy Controls in Agriculture. Agronomy 2021, 11, 110. [Google Scholar] [CrossRef]
  35. Chu, L.; Liu, Q.; Huang, C.; Liu, G. Monitoring of winter wheat distribution and phenological phases based on MODIS time-series: A case study in the Yellow River Delta, China. J. Integr. Agric. 2016, 15, 2403–2416. [Google Scholar] [CrossRef]
  36. Qiu, B.; Luo, Y.; Tang, Z.; Chen, C.; Lu, D.; Huang, H.; Chen, Y.; Chen, N.; Xu, W. Winter wheat mapping combining variations before and after estimated heading dates. ISPRS J. Photogramm. 2017, 123, 35–46. [Google Scholar] [CrossRef]
  37. Jin, M.; Zhang, R.; Sun, L.; Gao, Y. Temporal and spatial soil water management: A case study in the Heilonggang region, PR China. Agric. Water Manag. 1999, 42, 173–187. [Google Scholar] [CrossRef]
  38. He, J.; Li, H.; Rasaily, R.G.; Wang, Q.; Cai, G.; Su, Y.; Qiao, X.; Liu, L. Soil properties and crop yields after 11 years of no tillage farming in wheat–maize cropping system in North China Plain. Soil Tillage Res. 2011, 113, 48–54. [Google Scholar] [CrossRef]
  39. Jia, K.; Wu, B.; Li, Q. Crop classification using HJ satellite multispectral data in the North China Plain. J. Appl. Remote Sens. 2013, 3576, 073576. [Google Scholar] [CrossRef]
  40. Jingsong, S.; Guangsheng, Z.; Xinghua, S. Climatic suitability of the distribution of the winter wheat cultivation zone in China. Eur. J. Agron. 2012, 43, 77–86. [Google Scholar] [CrossRef]
  41. Kang, S.; Eltahir, E.A.B. North China Plain threatened by deadly heatwaves due to climate change and irrigation. Nat Commun 2018, 9, 2894. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  42. Wang, X.; Li, X.; Xin, L. Impact of the shrinking winter wheat sown area on agricultural water consumption in the Hebei Plain. J. Geogr. Sci. 2014, 24, 313–330. [Google Scholar] [CrossRef]
  43. Wu, Z.; Thenkabail, P.; Mueller, R.; Zakzeski, A.; Melton, F.; Johnson, L.; Rosevelt, C.; Dwyer, J.; Jones, J.; Verdin, J. Seasonal cultivated and fallow cropland mapping using MODIS-based automated cropland classification algorithm. J. Appl. Remote Sens. 2014, 8, 083685. [Google Scholar] [CrossRef] [Green Version]
  44. Chen, Y.; Lu, D.; Moran, E.; Batistella, M.; Dutra, L.V.; Sanches, I.D.A.; da Silva, R.F.B.; Huang, J.; Luiz, A.J.B.; de Oliveira, M.A.F. Mapping croplands, cropping patterns, and crop types using MODIS time-series data. Int. J. Appl. Earth Obs. 2018, 69, 133–147. [Google Scholar] [CrossRef]
  45. Wardlow, B.D.; Egbert, S.L. A comparison of MODIS 250-m EVI and NDVI data for crop mapping: A case study for southwest Kansas. Int. J. Remote Sens. 2010, 31, 805–830. [Google Scholar] [CrossRef]
  46. Onojeghuo, A.; Blackburn, G.; Wang, Q.; Atkinson, P.; Kindred, D.; Miao, Y. Rice crop phenology mapping at high spatial and temporal resolution using downscaled MODIS time-series. GISci. Remote Sens. 2018, 55, 659–677. [Google Scholar] [CrossRef] [Green Version]
  47. Yang, J.; Huang, X. The 30 m annual land cover dataset and its dynamics in China from 1990 to 2019. Earth Syst. Sci. Data 2021, 13, 3907–3925. [Google Scholar] [CrossRef]
  48. Savitzky, A.; Golay, M.J.E. Smoothing and Differentiation of Data by Simplified Least Squares Procedures. Anal. Chem. 1964, 36, 1627–1639. [Google Scholar] [CrossRef]
  49. Chen, J.; Jönsson, P.; Tamura, M.; Gu, Z.; Matsushita, B.; Eklundh, L. A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky–Golay filter. Remote Sens. Environ. 2004, 91, 332–344. [Google Scholar] [CrossRef]
  50. Wardlow, B.D.; Egbert, S.L. Large-area crop mapping using time-series MODIS 250 m NDVI data: An assessment for the U.S. Central Great Plains. Remote Sens. Environ. 2008, 112, 1096–1116. [Google Scholar] [CrossRef]
  51. Pan, Z.; Huang, J.; Zhou, Q.; Wang, L.; Cheng, Y.; Zhang, H.; Blackburn, G.A.; Yan, J.; Liu, J. Mapping crop phenology using NDVI time-series derived from HJ-1 A/B data. Int. J. Appl. Earth Obs. 2015, 34, 188–197. [Google Scholar] [CrossRef] [Green Version]
  52. Jonsson, P.; Eklundh, L. Seasonality extraction by function fitting to time-series of satellite sensor data. IEEE Trans. Geosci. Remote. Sens. 2002, 40, 1824–1832. [Google Scholar] [CrossRef]
  53. Leprieur, C.; Verstraete, M.M.; Pinty, B. Evaluation of the performance of various vegetation indices to retrieve vegetation cover from AVHRR data. Remote Sens. Rev. 1994, 10, 265–284. [Google Scholar] [CrossRef]
  54. Dong, C.; Zhao, G.; Qin, Y.; Wan, H. Area extraction and spatiotemporal characteristics of winter wheat-summer maize in Shandong Province using NDVI time series. PLoS ONE 2019, 14, e0226508. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  55. Yang, Y.; Wu, T.; Wang, S.; Li, J.; Muhanmmad, F. The NDVI-CV Method for Mapping Evergreen Trees in Complex Urban Areas Using Reconstructed Landsat 8 Time-Series Data. Forests 2019, 10, 139. [Google Scholar] [CrossRef] [Green Version]
  56. Weiss, E.; Marsh, S.; Pfirman, E. Application of NOAA-AVHRR NDVI time-series data to assess changes in Saudi Arabia’s rangelands. Int. J. Remote. Sens. 2001, 22, 1005–1027. [Google Scholar] [CrossRef]
  57. Fan, J.; Wu, B. A methodology for retrieving cropping index from NDVI profile. J. Remote Sens. 2004, 8, 628–636. [Google Scholar]
  58. Maxwell, S.K.; Sylvester, K.M. Identification of “ever-cropped” land (1984–2010) using Landsat annual maximum NDVI image composites: Southwestern Kansas case study. Remote Sens. Environ. 2012, 121, 186–195. [Google Scholar] [CrossRef] [Green Version]
  59. Goldblatt, R.; Stuhlmacher, M.F.; Tellman, B.; Clinton, N.; Hanson, G.; Georgescu, M.; Wang, C.; Serrano-Candela, F.; Khandelwal, A.K.; Cheng, W.-H.; et al. Using Landsat and nighttime lights for supervised pixel-based image classification of urban land cover. Remote Sens. Environ. 2018, 205, 253–275. [Google Scholar] [CrossRef]
  60. Janssen, L.L.F.; Vanderwel, F.J.M. Accuracy assessment of satellite derived land-cover data: A review. Photogramm. Eng. Remote Sens. 1994, 60, 419–426. [Google Scholar]
  61. Lu, S.; Deng, R.; Liang, Y.; Xiong, L.; Ai, X.; Qin, Y. Remote Sensing Retrieval of Total Phosphorus in the Pearl River Channels Based on the GF-1 Remote Sensing Data. Remote Sens. 2020, 12, 1420. [Google Scholar] [CrossRef]
  62. Akinyemi, F.O.; Pontius, R.G.; Braimoh, A.K. Land change dynamics: Insights from Intensity Analysis applied to an African emerging city. J. Spat. Sci. 2017, 62, 69–83. [Google Scholar] [CrossRef]
  63. Foody, G.M. Explaining the unsuitability of the kappa coefficient in the assessment and comparison of the accuracy of thematic maps obtained by image classification. Remote Sens. Environ. 2020, 239, 111630. [Google Scholar] [CrossRef]
  64. Pontius, R.G. Component intensities to relate difference by category with difference overall. Int. J. Appl. Earth Obs. 2019, 77, 94–99. [Google Scholar] [CrossRef]
  65. Pontius, R.G.; Millones, M. Death to Kappa: Birth of quantity disagreement and allocation disagreement for accuracy assessment. Int. J. Remote Sens. 2011, 32, 4407–4429. [Google Scholar] [CrossRef]
  66. Singh, A. Review Article Digital change detection techniques using remotely-sensed data. Int. J. Remote Sens. 1989, 10, 989–1003. [Google Scholar] [CrossRef] [Green Version]
  67. Liu, Y.; Lu, S.; Lu, X.; Wang, Z.; Chen, C.; He, H. Classification of Urban Hyperspectral Remote Sensing Imagery Based on Optimized Spectral Angle Mapping. J. Indian Soc. Remote. 2019, 47, 289–294. [Google Scholar] [CrossRef]
  68. Melesse, A.; Jordan, J. A comparison of fuzzy vs. augmented-ISODATA classification algorithms for cloud-shadow discrimination from Landsat images. Photogramm. Eng. Remote Sens. 2002, 68, 905–911. [Google Scholar] [CrossRef]
  69. Qiu, S.; Zhu, Z.; He, B. Fmask 4.0: Improved cloud and cloud shadow detection in Landsats 4–8 and Sentinel-2 imagery. Remote Sens. Environ. 2019, 231, 111205. [Google Scholar] [CrossRef]
  70. Hunter, E.L.; Power, C.H. An assessment of two classification methods for mapping Thames Estuary intertidal habitats using CASI data. Int. J. Remote Sens. 2002, 23, 2989–3008. [Google Scholar] [CrossRef]
  71. Shahriar Pervez, M.; Budde, M.; Rowland, J. Mapping irrigated areas in Afghanistan over the past decade using MODIS NDVI. Remote Sens. Environ. 2014, 149, 155–165. [Google Scholar] [CrossRef] [Green Version]
  72. Ahmad, A.; Sufahani, S.F. Analysis of Landsat 5 TM data of Malaysian land covers using ISODATA clustering technique. In Proceedings of the 2012 IEEE Asia-Pacific Conference on Applied Electromagnetics (APACE), Melaka, Malaysia, 11–13 December 2012; pp. 92–97. [Google Scholar] [CrossRef] [Green Version]
  73. Li, J.; Lei, H. Tracking the spatio-temporal change of planting area of winter wheat-summer maize cropping system in the North China Plain during 2001–2018. Comput. Electron. Agric. 2021, 187, 106222. [Google Scholar] [CrossRef]
  74. Xie, H.; Cheng, L.; Lv, T. Factors Influencing Farmer Willingness to Fallow Winter Wheat and Ecological Compensation Standards in a Groundwater Funnel Area in Hengshui, Hebei Province, China. Sustainability 2017, 9, 839. [Google Scholar] [CrossRef] [Green Version]
  75. Ju, J.; Masek, J.G. The vegetation greenness trend in Canada and US Alaska from 1984–2012 Landsat data. Remote Sens. Environ. 2016, 176, 1–16. [Google Scholar] [CrossRef]
  76. Sakamoto, T.; Yokozawa, M.; Toritani, H.; Shibayama, M.; Ishitsuka, N.; Ohno, H. A crop phenology detection method using time-series MODIS data. Remote Sens. Environ. 2005, 96, 366–374. [Google Scholar] [CrossRef]
  77. Li, F.; Ren, J.; Wu, S.; Zhao, H.; Zhang, N. Comparison of Regional Winter Wheat Mapping Results from Different Similarity Measurement Indicators of NDVI Time Series and Their Optimized Thresholds. Remote Sens. 2021, 13, 1162. [Google Scholar] [CrossRef]
  78. Zhang, X.; Zhang, M.; Zheng, Y.; Wu, B. Crop Mapping Using PROBA-V Time Series Data at the Yucheng and Hongxing Farm in China. Remote Sens. 2016, 8, 915. [Google Scholar] [CrossRef] [Green Version]
  79. Junges, A.; Fontana, D.; Pinto, D. Identification of croplands of winter cereals in Rio Grande do Sul state, Brazil, through unsupervised classification of normalized difference vegetation index images. Eng. Agrícola 2013, 33, 883–895. [Google Scholar] [CrossRef] [Green Version]
  80. Kar, S.; Rathore, V.S.; Champati Ray, P.K.; Sharma, R.; Swain, S.K. Classification of river water pollution using Hyperion data. J. Hydrol. 2016, 537, 221–233. [Google Scholar] [CrossRef]
  81. Nguyen, L.H.; Joshi, D.R.; Clay, D.E.; Henebry, G.M. Characterizing land cover/land use from multiple years of Landsat and MODIS time series: A novel approach using land surface phenology modeling and random forest classifier. Remote Sens. Environ. 2020, 238, 111017. [Google Scholar] [CrossRef]
  82. Wang, L.; Qi, F.; Shen, X.; Huang, J. Monitoring Multiple Cropping Index of Henan Province, China Based on MODIS-EVI Time Series Data and Savitzky-Golay Filtering Algorithm. Comput. Modeling Eng. Sci. 2019, 119, 331–348. [Google Scholar] [CrossRef] [Green Version]
  83. Yang, J.; Guo, N.; Huang, L.; Jia, J. Ananlyses on MODIS-NDVI Index Saturation in Northwest China. Plateau Meteorol. 2008, 27, 896–903. [Google Scholar]
  84. Faramarzi, M.; Heidarizadi, Z.; Mohamadi, A.; Heydari, M. Detection of vegetation cover changes using normalized difference vegetation index in semi-arid rangeland in western Iran. J. Agric. Sci. Technol. 2017, 20, 51–60. [Google Scholar]
  85. Chaves, M.E.D.; Picoli, M.C.A.; Sanches, I.D. Recent Applications of Landsat 8/OLI and Sentinel-2/MSI for Land Use and Land Cover Mapping: A Systematic Review. Remote Sens. 2020, 12, 3062. [Google Scholar] [CrossRef]
  86. Pastick, N.J.; Wylie, B.K.; Wu, Z. Spatiotemporal Analysis of Landsat-8 and Sentinel-2 Data to Support Monitoring of Dryland Ecosystems. Remote Sens. 2018, 10, 791. [Google Scholar] [CrossRef] [Green Version]
  87. Zeng, L.; Wardlow, B.D.; Xiang, D.; Hu, S.; Li, D. A review of vegetation phenological metrics extraction using time-series, multispectral satellite data. Remote Sens. Environ. 2020, 237, 111511. [Google Scholar] [CrossRef]
  88. Pastick, N.J.; Dahal, D.; Wylie, B.K.; Parajuli, S.; Boyte, S.P.; Wu, Z. Characterizing Land Surface Phenology and Exotic Annual Grasses in Dryland Ecosystems Using Landsat and Sentinel-2 Data in Harmony. Remote Sens. 2020, 12, 725. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Location of the study area in (a) China and (b) Hebei Province. (c) The land-use types obtained from the annual China Land Cover Dataset (CLCD) in 2017.
Figure 1. Location of the study area in (a) China and (b) Hebei Province. (c) The land-use types obtained from the annual China Land Cover Dataset (CLCD) in 2017.
Remotesensing 13 04870 g001
Figure 2. Main phenological periods of different crops in the study area.
Figure 2. Main phenological periods of different crops in the study area.
Remotesensing 13 04870 g002
Figure 3. Flowchart of winter wheat mapping and accuracy assessment. NDVI, normalized difference vegetation index; S-G, Savitzky-Golay; NDVIσ, standard deviation of NDVI; NDVIμ, mean value of NDVI; NDVI_COVfp, variation coefficient of time-series NDVI in fluctuation period; PI, peak index; SDI, slope difference index; SAM, Spectral Angle Mapping; ISOData, Spectral Angle Mapping; CBAH, Combining variations Before and After estimated Heading dates; WWI, Winter Wheat Index.
Figure 3. Flowchart of winter wheat mapping and accuracy assessment. NDVI, normalized difference vegetation index; S-G, Savitzky-Golay; NDVIσ, standard deviation of NDVI; NDVIμ, mean value of NDVI; NDVI_COVfp, variation coefficient of time-series NDVI in fluctuation period; PI, peak index; SDI, slope difference index; SAM, Spectral Angle Mapping; ISOData, Spectral Angle Mapping; CBAH, Combining variations Before and After estimated Heading dates; WWI, Winter Wheat Index.
Remotesensing 13 04870 g003
Figure 4. Research area map and sample data from Landsat-8 OLI images in Heilonggang region. The ranges of the sample areas are indicated by yellow lines. (iiii) represent the OLI images at overwintering, heading and harvesting stages of winter wheat, respectively. Rectangle a is characterized by build-up land; rectangle b is characterized by winter wheat; rectangle c is characterized by other vegetation.
Figure 4. Research area map and sample data from Landsat-8 OLI images in Heilonggang region. The ranges of the sample areas are indicated by yellow lines. (iiii) represent the OLI images at overwintering, heading and harvesting stages of winter wheat, respectively. Rectangle a is characterized by build-up land; rectangle b is characterized by winter wheat; rectangle c is characterized by other vegetation.
Remotesensing 13 04870 g004
Figure 5. (a) NDVI time series profiles before and after Savitzky-Golay (S-G) filter; (b) the date of MODIS images corresponding to each phase. The letters below the time-series profiles in (a) denote different phenological phases of winter wheat. A: overwintering, B: greening to jointing, C: heading to grouting, D: maturing to harvesting. DOY, the day of a year.
Figure 5. (a) NDVI time series profiles before and after Savitzky-Golay (S-G) filter; (b) the date of MODIS images corresponding to each phase. The letters below the time-series profiles in (a) denote different phenological phases of winter wheat. A: overwintering, B: greening to jointing, C: heading to grouting, D: maturing to harvesting. DOY, the day of a year.
Remotesensing 13 04870 g005
Figure 6. The mean NDVI time series profiles of winter wheat and other land cover types in 2017. PI1 and PI2 is characterized the peak location of NDVI time-series curve of winter wheat; SDI is characterized the slope of NDVI time-series curve of winter wheat in maturity-harvest stage; NDVI_COVfp is characterized by the variation coefficient of time-series NDVI in fluctuation period. DOY, the day of a year.
Figure 6. The mean NDVI time series profiles of winter wheat and other land cover types in 2017. PI1 and PI2 is characterized the peak location of NDVI time-series curve of winter wheat; SDI is characterized the slope of NDVI time-series curve of winter wheat in maturity-harvest stage; NDVI_COVfp is characterized by the variation coefficient of time-series NDVI in fluctuation period. DOY, the day of a year.
Remotesensing 13 04870 g006
Figure 7. The coefficient of variation (NDVI_COVfp) value and slope difference index (SDI) value corresponding to each curve (right).
Figure 7. The coefficient of variation (NDVI_COVfp) value and slope difference index (SDI) value corresponding to each curve (right).
Remotesensing 13 04870 g007
Figure 8. Model threshold analysis results: variation of (a) NDVI_COVfp, (b) SDI, (c) CBAH, and (d) WWI with changes in proportion of winter wheat in mixed pixels.
Figure 8. Model threshold analysis results: variation of (a) NDVI_COVfp, (b) SDI, (c) CBAH, and (d) WWI with changes in proportion of winter wheat in mixed pixels.
Remotesensing 13 04870 g008
Figure 9. Spatial distribution of winter wheat using COV_PSDI method in and around Heilonggang from 2014 to 2017 (Reference system is WGS84/UTM zone 49N).
Figure 9. Spatial distribution of winter wheat using COV_PSDI method in and around Heilonggang from 2014 to 2017 (Reference system is WGS84/UTM zone 49N).
Remotesensing 13 04870 g009
Figure 10. Comparison of different winter wheat mapping methods. (a) 30-m Landsat-8 OLI image combined as false color (5/4/3) on 8 June 2017. Extraction results using (b) COV_PSDI, (c) SAM, (d) ISODATA, (e) CBAH, and (f) WWI.
Figure 10. Comparison of different winter wheat mapping methods. (a) 30-m Landsat-8 OLI image combined as false color (5/4/3) on 8 June 2017. Extraction results using (b) COV_PSDI, (c) SAM, (d) ISODATA, (e) CBAH, and (f) WWI.
Remotesensing 13 04870 g010
Figure 11. Annual variation coefficient of NDVI time-series profile (NDVI_COVannual) corresponding to each land cover type.
Figure 11. Annual variation coefficient of NDVI time-series profile (NDVI_COVannual) corresponding to each land cover type.
Remotesensing 13 04870 g011
Table 1. Information and Sources of Remote Sensing and Land Cover Data.
Table 1. Information and Sources of Remote Sensing and Land Cover Data.
Satellite DataAcquisition TimePath/RowSpatial ResolutionSpectral ResolutionProduct GradeSource
MODIS 1 (384 images)1 November 2013–19 December 2017h26v04
h26v05
h27v04
h27v05
250 mNDVIMOD13Q1https://ladsweb.modaps.eosdis.nasa.gov/ (accessed on 25 April 2020)
Landsat-8 (12 images)15 June 2017124/03330 mGreen: (μm)
(0.525–0.600)
OLI 3
Reflectance
http://earthexplorer.usgs.gov/ (accessed on 10 April 2020)
14 May 2017
5 December 2016
15 June 2017124/034Red: (μm)
(0.630–0.680)
28 April 2017
22 January 2017
8 June 2017123/034
123/035
NIR 2: (μm)
(0.845–0.885)
7 May 2017
14 December 2016
CLCD 4 [47]2017 30 mLand CoverVersion 1.0.0https://doi.org/10.5281/zenodo.4417810 (accessed on 10 March 2021)
1 MODIS, Moderate Resolution Imaging Spectroradiometer; 2 NIR, near-infrared; 3 OLI, Operational Land Imager; 4 CLCD, annual China Land Cover Dataset.
Table 2. Parameters of WWI and CBAH models.
Table 2. Parameters of WWI and CBAH models.
ModelEquationReference
CBAH 1 { E V E 2 = E V I 2 m a x 1 E V I 2 m i n 1 E V L 3 = E V I 2 m a x 1 E V I 2 m i n 2
EVI2max1 is the maximum EVI value in the heading stage; EVI2min1 and EVI2min2 is the minimum EVI value in the seeding and harvesting stage, respectively.
Qiu et al. [36]
WWI 4 W W I = { ( N D V I m a x 1 N D V I m i n 1 ) × ( N D V I m a x 2 N D V I m i n 2 ) ,   i f   N D V I m a x 1 < N D V I m i n 1   a n d   N D V I m a x 2 < N D V I m i n 2 o r   N D V I m i n 1 < N D V I m a x 1 < 0   a n d   N D V I m i n 2 < N D V I m a x 2 < 0 ( N D V I m a x 1 N D V I m i n 1 ) × ( N D V I m a x 2 N D V I m i n 2 ) ,   o t h e r w i s e .
NDVImax1 and NDVImax2 are the maximum NDVI value in the over-wintering stage and heading stage, respectively; NDVImin1 and NDVImin2 are the minimum NDVI value in the seeding stage and harvesting stage, respectively.
Qu et al. [33]
1 CBAH, Combining variations before and after estimated heading dates; 2 EVE, early growth stage (seeding stage to heading stage); 3 EVL, late growth stage (heading stage to harvesting stage); 4 WWI, winter wheat index.
Table 3. Comparison of the mapping accuracy of different models in 2017.
Table 3. Comparison of the mapping accuracy of different models in 2017.
ModelsIdentificationVerificationUA
Winter WheatNon-Winter Wheat
COV_PSDIWinter wheat4643193.74%
Non-winter wheat28477
PA94.32%
OA94.10%
SAMWinter wheat4515289.66%
Non-winter wheat41456
PA91.67%
OA90.70%
ISODataWinter wheat4578883.85%
Non-winter wheat35420
PA92.89%
OA87.70%
CBAHWinter wheat4616188.31%
Non-winter wheat31447
PA93.70%
OA90.80%
WWIWinter wheat48910782.05%
Non-winter wheat3401
PA99.39%
OA89.00%
Table 4. PE for different classification methods of 2014–2017 (%).
Table 4. PE for different classification methods of 2014–2017 (%).
YearStatistical Data (km2)COV_PSDISAMISODataCBAHWWI
Extracted Area (km2)PE (%)Extracted Area (km2)PE (%)Extracted Area (km2)PE (%)Extracted Area (km2)PE (%)Extracted Area (km2)PE (%)
201421,615.022,076.22.120,776.73.927,375.326.620,961.33.023,123.07.0
201521,297.222,072.83.622,718.76.729,249.937.315,938.425.224,525.615.2
201620,739.720,416.81.622,272.77.421,640.64.316,243.921.724,954.620.3
201722,062.822,821.33.424,469.810.926,944.922.124,681.911.923,875.38.2
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Zhang, X.; Liu, K.; Wang, S.; Long, X.; Li, X. A Rapid Model (COV_PSDI) for Winter Wheat Mapping in Fallow Rotation Area Using MODIS NDVI Time-Series Satellite Observations: The Case of the Heilonggang Region. Remote Sens. 2021, 13, 4870. https://doi.org/10.3390/rs13234870

AMA Style

Zhang X, Liu K, Wang S, Long X, Li X. A Rapid Model (COV_PSDI) for Winter Wheat Mapping in Fallow Rotation Area Using MODIS NDVI Time-Series Satellite Observations: The Case of the Heilonggang Region. Remote Sensing. 2021; 13(23):4870. https://doi.org/10.3390/rs13234870

Chicago/Turabian Style

Zhang, Xiaoyuan, Kai Liu, Shudong Wang, Xin Long, and Xueke Li. 2021. "A Rapid Model (COV_PSDI) for Winter Wheat Mapping in Fallow Rotation Area Using MODIS NDVI Time-Series Satellite Observations: The Case of the Heilonggang Region" Remote Sensing 13, no. 23: 4870. https://doi.org/10.3390/rs13234870

APA Style

Zhang, X., Liu, K., Wang, S., Long, X., & Li, X. (2021). A Rapid Model (COV_PSDI) for Winter Wheat Mapping in Fallow Rotation Area Using MODIS NDVI Time-Series Satellite Observations: The Case of the Heilonggang Region. Remote Sensing, 13(23), 4870. https://doi.org/10.3390/rs13234870

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