Maps of cropping patterns in China during 2015–2021

Study area

There is a long history of diversified cropping patterns due to the climatic and topographic complexity in China⁴. Cropping intensity increases from north to south, and multiple cropping dominates in regions south of 40⁰N⁴. For example, multiple cropping systems of double rice and winter wheat plus maize are popular in the Middle-lower Yangtze river plain and the Huang-Huai-Hai plain, respectively (Fig. 1)²². Three staple crops, maize, paddy rice, and wheat, are widely distributed across the country (Figure S1). These three major crops contributed to more than half (57.08%) of the total sown area in China in 2020 (http://www.stats.gov.cn/english/).

MODIS images and pre-processing

We used the 500 m 8-day composite Moderate Resolution Imaging Spectroradiometer (MODIS) surface reflectance products (MOD09A1) from 2015 to 2021. Three spectral indices were calculated: the 2-band Enhanced Vegetation Index (EVI2)²³, LSWI¹⁶, and Normalized Multi-band Drought Index (NMDI)²⁴ (Fig. 2). The functions of EVI2, LSWI, and NMDI are provided in Eqs. 1–3 as follows.

where, ρ_NIR, ρ_Red, ρ_SWIR6 and ρ_SWIR7 represented the surface reflectance values from the red (620–670 nm), Near-infrared (841–875 nm), short wave infrared band centered at 1640 nm (1628–1652 nm) and 2130 nm (2105–2155 nm), respectively.

For each spectral index (EVI2, LSWI, and NMDI), a daily continuous time series was developed based on the cloud-free observations using the Whittaker Smoother (WS)²⁵. The WS smoother performed well in multiple cropping regions and therefore was applied here²⁶.

Validation data and other datasets

The validation data in this study included the ground truth reference data and agricultural census data. The ground truth reference data were collected in major agricultural regions with GPS receivers and digital cameras during the study period (2015–2021) (Fig. 1, Table S1). For each sampling site, the geographic location and crop types were recorded. The reliability of ground survey data was improved through visual confirmation based on high-resolution images in Google Earth. Some reference sites with small field sizes were removed to considering the mixed-pixel problems of MODIS images. Finally, we obtained a total of 18,379 ground samples collected during 2015–2021 (Table S1). All the ground truth reference data were used to validate the cropping pattern data in its corresponding year. Agricultural census data were obtained from the National Statistical Bureau of China (NSBC) (http://www.stats.gov.cn/english/), which was collected through sampling statistics. The crop calendar data from agro-meteorological stations recorded the sowing, seedling, tillering, heading, and harvesting dates of winter wheat (210 sites) or spring wheat (90 sites). The calendar data were applied to establish the trend surfaces of key phenological stages of winter wheat and spring wheat, respectively. The crop calendar data were provided by the National Meteorological Information Center, China Meteorological Administration.

The cropland distribution data were derived from the 30 m GlobeLand30 global land cover data in 2020²⁷. The total accuracy of GlobeLand30 in 2020 is 85.72%, and the Kappa coefficient is 0.82 (www.globallandcover.com). To correspond to MODIS images, the 30 m cropland pixels from GlobeLand30 data were spatially aggregated to a 500 m cropland fraction map. For simplification, we classified pixel purity of MODIS pixels into three groups: cropland percentages of >90%, 50–90%, and <50% were labeled as pure, moderate mixed, and seriously mixed pixels, respectively. MODIS pixels with very small cropland fraction (i. e. <30%) were not accounted. The pure, moderate mixed and seriously mixed groups occupied 39%, 42%, and 19% of MODIS pixels in China, respectively (Figure S2).

Overview of the cropping pattern mapping approach

A cropping pattern is referred to as the yearly sequence and spatial distribution of major crops on a specific piece of cropland. Consequently, cropping pattern mapping should provide information on cropping intensity as well as crop types. When multiple cropping is cultivated, we derived the plantation sequence of two or triple crops. For example, the cropping pattern of “winter wheat-maize” represents double cropping with a rotation of winter wheat plus maize. We conducted cropping pattern mapping processes using MATLAB software (Fig. 2). Annual cropping pattern maps were obtained by deriving cropping intensity as well as mapping three staple crop types (paddy rice, maize, wheat). These knowledge-based mapping algorithms were described in the following sections (Fig. 2).

Deriving cropping intensity

The vegetation indices (VI) peak-based algorithms have been widely applied to identify cropping intensity in previous studies²⁸. However, the VI peak-based algorithms were challenged by the changes of VI temporal profiles in different cropping patterns over large areas and multiple years^15,17. An automatic Cropping Intensity extraction method based on the Isolines of Wavelet Spectra (CIIWS) was proposed with considerations of complex intra-class variability of VI temporal profiles^29,30. The cropping index was identified based on three main features, the skeleton width, maximum number of strong brightness centers, and the intersection of their scale intervals, derived from wavelet spectra (Fig. 2)^29,30. The CIIWS cropping intensity mapping algorithm is capable of deriving cropping intensity automatically, which is robust to intra-class variability such as the phenology shift, strengthened or lessened crop growth, or crop diversity^29,30. The wavelet features-based cropping intensity mapping algorithm was applied in mainland China to derive annual cropping intensity from 1982 to 2013, with an overall accuracy of 91.63%⁴. Therefore, this wavelet features-based cropping intensity mapping algorithm was applied in this study^29,30.

Mapping paddy rice, maize and wheat based on phenological indicators

Algorithms for mapping paddy rice, maize and wheat were developed in related references^15,31,32 and applied in this study. We made critical improvements over these proposed crop mapping algorithms primarily in order to cope with the mixed-pixel problem in MODIS images^{15,29,30,31,32}. First, pixel purity-based thresholds were applied in the decision rules (Table 1) to cope with the mixed pixel problem of MODIS images. The target crop from mixed pixels was expected to show lower values than those from pure pixels when the target crop was highlighted by larger values in our proposed indicators (i.e., maize). The pixel purity-based thresholds were determined based on the accuracy assessment with agricultural census data in 2018. The pixel purity-based thresholds were applied to a national scale and multiple years (2015–2021) without adaptations. Second, the derived maps were further improved by incorporating their corresponding suitable areas. Specifically, the suitable areas of single rice and spring maize were estimated based on topographic and climatic conditions (elevation and accumulated temperature greater than 10 degrees)³³. Finally, we calculated the estimated areas of three staple crops from MODIS-derived products using the cropland fraction map. A concise description on the original mapping algorithms were provided in the following sections.

A phenology-based rice mapping algorithm was proposed through Combined Consideration of Vegetation phenology and Surface water variations (CCVS)³¹. Variation of LSWI in rice fields was smaller than that in other crops fields during the period from tillering to heading dates³¹. Therefore, the Ratio of Change amplitude of LSWI to 2-band Enhanced Vegetation Index 2 (RCLE) during that period was utilized as the primary metric for paddy rice mapping (Fig. 2, Table 1). The CCVS rice mapping algorithm was successfully applied in 15 provincial-level administrative units of southern China, which obtained an overall accuracy of 93–94%³¹. The CCVS rice mapping algorithm proved to be robust in terms of intra-class rice variability¹⁴.

An automatic Maize mapping was recently proposed through Exploring Leaf moisture variation during flowering Stage (MELS)³². One unique indicator for maize mapping was the Ratio of Cumulative Positive slope to Negative slope (RCPN) of NMDI during the flowering stage (Fig. 2, Table 1). Maize sites were highlighted by consistently higher values in our proposed phenology-based indicator (RCPN) compared to other crops. A simple rule was applied to derive maize³² (Table 1). The capability of the MELS method was verified in mainland China, with an overall accuracy of 91%³².

A phenology-based winter wheat mapping algorithm through Combining variations Before and After estimated Heading dates (CBAH) (Fig. 2) was exploited in this study¹⁵. The CBAH algorithm demonstrated the adaptability of VI temporal profiles to intra-class variability³⁴. This mapping algorithm was exploited in North China plain from 2001 to 2013, with an overall accuracy of around 90%¹⁵. The CBAH algorithm was extended to map both winter wheat and spring wheat at the national scale. The potential distribution of winter wheat and spring wheat was determined based on latitudes and the accumulated temperature above 5 °C with references to crop calendar data from agro-meteorological stations. And then the trend surfaces of key phenological stages (heading dates and early growing length) were established for winter wheat and spring wheat, respectively³⁵. For winter wheat, altitude and latitude were applied; and for spring wheat, the accumulated temperature above 5 °C and latitude were applied. Two phenology-based indicators were developed by exploring the variations of VI during the estimated early and late growth stages (Fig. 2, Table 1). A simple decision rule could be applied to identify wheat by these two phenology-based indicators (Table 1).

Source: Ecology - nature.com