Multi-year monitoring of paddy rice planting area in Northeast China using MODIS time series data

Multi-year monitoring of paddy rice planting area in Northeast China using MODIS time series data 1 2 first-author 1 2 3 corresp ‡Corresponding Author ‡Corresponding Author †E-mail:hjf@zju.edu.cn †E-mail:hjf@zju.edu.cn hjf@zju.edu.cn 1 2 Institute of Agricultural Remote Sensing and Information Application, Zhejiang University, Hangzhou 310058, China Institute of Agricultural Remote Sensing and Information Application, Zhejiang University, Hangzhou 310058, China Key Laboratory of Agricultural Remote Sensing and Information System of Zhejiang Province, Hangzhou 310058, China Key Laboratory of Agricultural Remote Sensing and Information System of Zhejiang Province, Hangzhou 310058, China Ministry of Education Key Laboratory of Environmental Remediation and Ecological Health, Zhejiang University, Hangzhou 310058, China Ministry of Education Key Laboratory of Environmental Remediation and Ecological Health, Zhejiang University, Hangzhou 310058, China The objective of this study was to investigate the tempo-spatial distribution of paddy rice in Northeast China using moderate resolution imaging spectroradiometer (MODIS) data. We developed an algorithm for detection and estimation of the transplanting and flooding periods of paddy rice with a combination of enhanced vegetation index (EVI) and land surface water index with a central wavelength at 2 130 nm (LSWI 2130 ). In two intensive sites in Northeast China, fine resolution satellite imagery was used to validate the performance of the algorithm at pixel and 3×3 pixel window levels, respectively. The commission and omission errors in both of the intensive sites were approximately less than 20%. Based on the algorithm, annual distribution of paddy rice in Northeast China from 2001 to 2009 was mapped and analyzed. The results demonstrated that the MODIS-derived area was highly correlated with published agricultural statistical data with a coefficient of determination ( R 2 ) value of 0.847. It also revealed a sharp decline in 2003, especially in the Sanjiang Plain located in the northeast of Heilongjiang Province, due to the oversupply and price decline of rice in 2002. These results suggest that the approaches are available for accurate and reliable monitoring of rice cultivated areas and variation on a large scale. Paddy rice Moderate resolution imaging spectroradiometer (MODIS) Northeast China Enhanced vegetation index Land surface water index Introduction 1 Rice ( Oryza sativa L.) is the primary staple food source for more than half of the world’s population and has profound influence on the livelihood of farmers (Khush, 2005 2005 ). Timely and reliable information about paddy rice distribution over wide areas is a topic of global interest and has special significance to China, where grain accounts for more than 40% of crop production (Shao et al., 2001 2001 ). Northeast China is the major japonica rice producing region in China with a grain yield accounting for 44.6% of the total output of national japonica rice in 2009 (Liu et al., 2010 2010 ). Mapping rice plant areas will also lead to more accurate assessments of irrigation water supply (Rosenzweig et al., 2004 2004 ), grain production (Amano et al., 1993 1993 ; Doraiswamy et al., 2005 2005 ), and net primary production (Ma, 2008 2008 ). In addition, the rice field has been recognized as one of the major anthropogenic emission sources of greenhouse gases (GHGs) and contributes to over 10% of the global atmospheric input of methane (Neue, 1993 1993 ). Monitoring and mapping of paddy rice is thus very important for food security assessment and planning, environmental sustainability, and government decision making (Sakamoto et al., 2005 2005 ; Zhang et al., 2009 2009 ; Gumma et al., 2011 2011 ). Steady progress has been made toward monitoring crop growth conditions and estimating crop yield with remote sensing techniques for several decades. Compared with the traditional census data acquired through a sample survey, remote sensing is more feasible and less labor-intensive for providing spatial distribution of paddy rice in real time with its large-area coverage and frequent revisits. In the 1980s, a large-area crop inventory experiment (LACIE) project using Landsat multispectral scanner (MSS) was a satisfactory and successful trial for crop survey with remote sensing data (Macdonald and Hall, 1980 1980 ). Okamoto and Fukuhara ( 1996 1996 ) developed a model to estimate the area of paddy rice fields in pixels from Landsat thematic mapper (TM) data and perform validation in Hokkaido and Miyagi Prefectures in Japan. Synthetic aperture radar (SAR) is anticipated to be the primary high-resolution remote sensing data for rice detection, especially in tropical and subtropical regions with long-term cloudy and rainy weather (Shao et al., 2001 2001 ; Zhang et al., 2009 2009 ). Satisfactory classification accuracy was obtained on a local scale with the aforementioned high spatial resolution satellite data. However, the temporal analysis of land cover was constrained by the lower temporal resolution, limited coverage extent, and high cost of the images (Peng et al., 2011 2011 ). In recent years, moderate resolution imaging spectroradiometer (MODIS) on board the Terra and Aqua satellites was widely applied to monitoring rice areas and detecting the phenology stage on a regional scale due to its advantages of moderate spatial resolution and high revisited periods (Sakamoto et al., 2005 2005 ; Xiao et al., 2005 2005 ; Sun et al., 2009 2009 ; Peng et al., 2011 2011 ). Paddy rice fields are usually prepared by flooding a few days before rice seedlings are transplanted and have at least one wet growing season in shallow floodwater (Neue, 1993 1993 ). Since the intermittent flooding in a paddy rice field could generate a distinct signature, identification of the flooding and transplanting periods is the foundation of rice planting area extraction (Malingreau, 1986 1986 ; Bachelet, 1995 1995 ; Xiao et al., 2006 2006 ). Spectral bands or spectral indices, which are sensitive to both water and vegetation, are necessary for monitoring changes in the mixture of surface water and green vegetation in paddy rice fields (Sari et al., 2010 2010 ). So far, a large number of spectral indices have been developed by combining two or more bands to enhance the contrast between target (e.g., vegetation, water) and background and reduce effects of atmosphere, environment, solar illumination geometry, and sensor viewing conditions. Normalized difference vegetation index (NDVI) is constructed based on the red and near infrared (NIR) bands, which are located in the strong chlorophyll absorption region and high reflectance plateau of vegetation canopies, respectively. It is correlated to the leaf area index and chlorophyll content and has been widely used in many ways, such as for estimating crop yield and net primary production, identifying crop types, and detecting land cover change (Quarmby et al., 1993 1993 ; Lunetta et al., 2006 2006 ). Enhanced vegetation index (EVI) was constructed to overcome the limitations of NDVI such as saturation in high biomass regions, and adjust residual atmospheric contamination and background reflectance (Huete et al., 1997 1997 ; 2002 2002 ). A potentially better way to estimate vegetation water content is to use indices based on infrared range (1 240–3 000 nm). The spectral index combining NIR with shortwave infrared (SWIR) bands has the capability of retrieving canopy water content (Ceccato et al., 2002a 2002a ; 2002b 2002b ). Gao ( 1996 1996 ) proposed a normalized difference water index (NDWI) for estimation of vegetation liquid water using remote sensing data. Chen et al. ( 2005 2005 ) tested the ability of estimating vegetation water contents of corn and soybean with the surface reflectance at 1 640 nm and 2 130 nm, namely NDWI 1640 and NDWI 2130 , respectively. Hereafter another nomenclature of NDWI named land surface water index (LSWI) centered at 1 640 nm was used for identifying water properties in the flooding and transplanting stages of paddy rice with the condition of LSWI+threshold ( T )≥NDVI or LSWI+ T ≥EVI (Xiao et al., 2005 2005 ; 2006 2006 ; Sun et al., 2009 2009 ). However, Kim et al. ( 2004 2004 ) reported that the water index LSWI using the SWIR band centered at 2 130 nm appeared to be more useful for detecting vegetation water status. The SWIR band centered at 2 130 nm is less affected by ozone and Rayleigh scattering, water vapor, and aerosols than that at 1 640 nm and has a good response to cumulative rainfall (Vermote et al., 1997 1997 ; Chandrasekar et al., 2010 2010 ). Therefore, the combined index LSWI using NIR and SWIR centered at 2 130 nm (viz. LSWI 2130 ) was used in this study. The objective of this study is to: (1) verify the potential of LSWI 2130 combined with EVI to identify the flooding and transplanting stages of paddy rice in Northeast China; (2) evaluate the performance of our algorithm in main rice producing regions of Northeast China with TM images; (3) analyze the tempo-spatial distribution of rice in Northeast China during the period 2001–2009. Materials and methods 1 Study area and data 2 Study area 3 Northeast China, encompassing three provinces of Heilongjiang, Jilin and Liaoning, is one of the most important agricultural production bases in China. The study area extends in longitude from 118°50′5″ E to 134°46′26″ E and in latitude from 38°43′15″ N to 53°33′39″ N with a territory of 7.9×10 5 km 2 of which nearly 30% is arable (Fig. 1 1 ). Annual precipitation in the study area varies widely from 480 mm in the west to 900 mm along the east coast. This area has a boreal climate with an annual temperature of around 1–7 °C and the frost-free period was about 80–180 d (Liu et al., 2005 2005 ). The unique pattern of vegetation distribution and different population diversities was formed in the whole region because of its specific geographical and meteorological features (Luo et al., 2009 2009 ; Zhang L.Z. et al., 2011 2011 ). Paddy rice is one of the staple crops in Northeast China with high annual yield and taste quality because it is cultivated as a single-season crop in the major alluvial plains with a gentle and flat topographic relief, abundant precipitation, and fertile soils without any organic matter amended (Zhang Y. et al., 2011 2011 ). Satellite data 3 The National Aeronautics and Space Administration (NASA) provides different MODIS products processed with standard correction algorithms to users for free (available at http://reverb.echo.nasa.gov/reverb/ http://reverb.echo.nasa.gov/reverb/ ). The MODIS/Terra surface reflectance 8-d product (MOD09A1) is composed so as to have the best observation for each pixel and the lowest value of band 3 (blue band) over the 8-d period and has already been systematically corrected for the effects of gases, thin cirrus clouds, and aerosols scattering and absorption (Vermote and Vermeulen, 1999 1999 ). MOD09A1 also includes two quality assessment (QA) datasets at pixel and band levels, which are essential for users to note the scientific quality of the products with regard to their intended performance and remove areas of persistent cloud (Roy et al., 2002 2002 ). Each tile covers an area of 1 200 km×1 200 km and is produced at a spatial resolution of approximately 500 m with sinusoidal projection. In this study, five tiles (h25v03, h26v03, h26v04, h27v04, h27v05) in the period 2001–2009 were necessary for covering the whole region of Northeast China. Mosaic and reprojection of downloaded MODIS data were preprocessed with MODIS reprojection tool (MRT) provided by the U.S. Geological Survey (USGS) Earth Resources Observation and Science (EROS) Center. In this study, we selected two intensive sites (A and B) as the test regions (Fig. 1 1 ). Cloud-free or nearly cloud-free Landsat TM data at 30-m resolution for bands 1, 2, 3, 4, 5, 7 plus 120 m for band 6 were obtained and used to evaluate the performance of our algorithm for extracting paddy rice with MODIS data and to validate their accuracy (Table 1 1 ). The downloaded TM images were resized and reprojected to Albers equal conic projection. A paddy rice binary map of these sites was obtained from TM images using the maximum likelihood method. Site No. Location Acquired date Path/row number A Heilongjiang June 14, 2007 114/28 B Liaoning June 17, 2007; July 19, 2007; Sept. 5, 2007 119/31 Ancillary data 3 Annual provincial agricultural census data of the rice sown area and producer price indices of rice were collected from the China Statistical Yearbook and China Rural Statistical Yearbook databases, respectively (available at http://tongji.cnki.net/ http://tongji.cnki.net/ ). Annual provincial agricultural census data were used as validation data to evaluate the performance of results derived from the MODIS algorithm. A digital administrative map of Northeast China was provided by the China National Fundamental Geographic Information System. Methodology for extracting paddy rice pixels 2 Spectral index calculation and reconstruction 3 With the intention of discriminating paddy rice planting pixels with MODIS data, four spectral indices, namely NDVI, EVI, LSWI 2130 , and normalized difference snow index (NDSI), were calculated using the following equations: $${\text{NDVI}} = ({\rho _{{\text{band2}}}} - {\rho _{{\text{band1}}}})/({{\rho _{{\text{band2}}}} + {\rho _{{\text{band1}}}}})$$ , $${\text{EVI}} = ({{\rho _{{\text{band2}}}} - {\rho _{{\text{band1}}}}})/({{\rho _{{\text{band2}}}} + {\text{6}}.0{\rho _{{\text{band1}}}} - {\text{7}}.{\text{5}}{\rho _{{\text{band3}}}} + {\text{1}}})$$ , $${\text{LSW}}{{\text{I}}_{{\text{2130}}}} = ({{\rho _{{\text{band2}}}} - {\rho _{{\text{band7}}}}})/({{\rho _{{\text{band2}}}} + {\rho _{{\text{band7}}}}})$$ , $${\text{NDSI}} = ({{\rho _{{\text{band4}}}} - {\rho _{{\text{band6}}}}})/({{\rho _{{\text{band4}}}} + {\rho _{{\text{band6}}}}})$$ , where ρ band1 , ρ band2 , ρ band3 , ρ band4 , ρ band6 , and ρ band7 are the surface reflectances of band 1 (red), band 2 (NIR), band 3 (blue), band 4 (green), band 6 (SWIR: 1 640 nm), and band 7 (SWIR: 2 130 nm), respectively. The potential of LSWI to estimate vegetation water content and rainfall accumulation makes it possible to identify the flooding and transplanting periods of paddy rice, of which the spectral characteristic is predominated by water. Cloud-contaminated pixels were extracted according to the QA description of MOD09A1. In addition, the pixels whose reflectance of the blue band was greater than 0.1, but not labeled as clouds in the cloud quality flag of MOD09A1 dataset, were removed as abnormal data because of the effect of thick clouds (Sakamoto et al., 2005 2005 ). NDSI was calculated with the green and SWIR (1 640 nm) bands (Eq. ( 4 4 )). Snow cover pixels were identified with the threshold NDSI≥0.40 and band 2 reflectance > 0.11 (Hall et al., 2002 2002 ). In order to avoid underestimating the result caused by missing data contaminated by cloud, we interpolated the missing data in NDVI and EVI time-series data with local maximum fitting (LMF): $${\text{VI}}{\prime _t} = {\text{min}}[{\text{max}}\left( {{\text{V}}{{\text{I}}_t}_{ - {\text{3}}},{\text{ V}}{{\text{I}}_t}_{ - {\text{2}}},{\text{ V}}{{\text{I}}_t}_{ - {\text{1}}},{\text{ V}}{{\text{I}}_t}} \right),{\text{max}}\left( {{\text{V}}{{\text{I}}_t}_{ + {\text{3}}},{\text{ V}}{{\text{I}}_t}_{ + {\text{2}}},{\text{ V}}{{\text{I}}_t}_{ + {\text{1}}},{\text{ V}}{{\text{I}}_t}} \right)]$$ , where VI′ t and VI t are the reconstructed and observed vegetation indices at time t , respectively (Wada and Ohira, 2004 2004 ). Characteristics of spectral indices in rice flooding and transplanting periods 3 Paddy rice in Northeast China is well irrigated and planted on flooded soil. During the flooding and transplanting periods, the spectral reflectance is mainly predominated by water since the field is a mixture of water of 2 to 15 cm depth and rice seedlings (Xiao et al., 2002 2002 ). In the test rice field, located in Liaozhong County, Liaoyang City, Liaoning Province (122°54′19″ E, 41°36′ N), the spectral indices such as LSWI 2130 and EVI were analyzed. According to the phenological data observed at the meteorological station, the transplanting period and heading stage of paddy rice in Liaozhong County began on May 17 and Aug. 1. As shown in Fig. 2 2 , the reconstructed seasonal EVI time profile had a single peak in Aug. 2, which was well consistent with the agricultural cropping system of single rice. LSWI 2130 varied significantly in water content and increased sharply a few days after the transplanting period. Henceforward, LSWI 2130 increased more slightly than EVI after about 20 d of the transplanting period with the subsequent rapid vegetation growth of the rice. The seasonal profile of LSWI 2130 kept increasing since the beginning of the flooding and transplanting periods and reached a maximum in early August. As a result of snow cover and melting, LSWI 2130 presented high value in early March and late December of the year. EVI was found to be higher than 0.5LSWI 2130 in the flooding and transplanting stages which lasted for about two to three 8-d (Fig. 2 2 ). Accordingly, we developed a hypothesis for identifying flooding and transplanting pixels: LSWI 2130 > 2EVI in the test rice field. Afterwards, a second procedure to determine the rice pixel was implemented using the assumption that the mean EVI of the 7th to 11th 8-d composites is greater than or equal to 0.45. In fact, since individual farmers have different planting schedules, the cropping pattern of paddy rice in Northeast China varied both spatially and temporally. According to the records of agricultural meteorological stations, the transplanting data of paddy rice in Northeast China began from May to June, so that for each MOD09A1 tile from May to June, the initial rice map was generated and then merged into one rice map. Supplementary mask layer for rice extraction algorithm 3 In addition, it is necessary to separate the persistent watershed, wetland, and forest from the extracted rice field. The persistent water pixels with hydrophytes can be misidentified as paddy rice using the above-mentioned algorithm. In order to exclude the overestimation of paddy rice caused by persistent water, we analyzed the temporal profiles of NDVI and LSWI 2130 , and assumed that a pixel was covered by water if NDVI < 0.1 and LSWI 2130 > 0.2. The pixel was identified as persistent water if it was found to be covered by water in twenty or more 8-d composites all through the year. Then the annual persistent water maps were generated and excluded from the paddy rice results. Furthermore, the flooded forest, shrub, and wetland should be separated from the paddy rice. The seasonal NDVI profile is deemed as an indicator of greenness of vegetation and has been widely used to define the key phenological phases of vegetation dynamics by identifying its starting, ending, and inflexion points (Zhang et al., 2003 2003 ). Forest and natural wetland vegetations in Northeast China have a long growing season period from mid-late April to mid October and a consistently high NDVI (or LSWI 2130 ) throughout the year. Moreover, the high NDVI (or LSWI 2130 ) of paddy rice merely lasts for a few 8-d. According to the comparison of NDVI and LSWI 2130 , we found that forest and wetland could be discriminated from paddy rice if the MODIS pixel satisfied the relation NDVI > 0.4 and NDVI−LSWI 2130 > 0.05 over 15 8-d composites during the year. Companion Landsat TM data processing and accuracy assessment 2 In order to evaluate the performance of our MODIS-based algorithm, the classification result of the companion TM imagery at each intensive site was degraded to the same spatial resolution as MODIS data using the pixel aggregate method and used as the rice reference map for further validation (Barker and Burelhach, 1992 1992 ; Fang and Liang, 2005 2005 ). The accuracy was used to describe the agreement of results with the rice reference map. An error confusion matrix was applied in assessing commission error (CE; percentage of falsely identified rice pixels in the MODIS results) and omission error (OE; percentage of true rice pixels not found in the rice reference map). The formulae used for accuracy assessment were determined as $${\text{CE}} = {N_{{\text{commit}}}}/{N_{{\text{MODIS}}}} \times 100\%$$ , $${\text{User'}} {\text{s }} {\text{accuracy}} = 100\% - {\text{CE}}$$ , $${\text{OE}} = {N_{{\text{omit}}}}/{N_{{\text{TM}}}} \times 100\%$$ , $${\text{Producer'}} {\text{s }} {\text{accuracy}} = 100\% - {\text{OE}}$$ , where N commit and N omit represent the amount of commit and omit paddy rice pixels, respectively, and N MODIS and N TM represent the amount of paddy rice pixels in the MODIS algorithm and aggregated rice reference map, respectively. The accuracy assessment was first applied at the pixel level, which means that each pixel in the MODIS-derived result was compared to the rice reference binary map. Furthermore, in order to compensate for potential edge effects caused by the spatial aggregation of the TM image and the geometrical mismatch between TM and MODIS data, the assessment was implement using the 3×3 pixel moving window (Zhan et al., 2000 2000 ; Fang and Liang, 2005 2005 ). Results 1 Accuracy assessment of paddy rice extraction algorithm with rice reference map in the intensive sites 2 At pixel level, the accuracy assessment was conducted per pixel with the comparison between MODIS-derived rice distribution and the rice reference map. In site A, the CE and OE were 26% and 49%, respectively, i.e., the user’s and producer’s accuracies were 74% and 51%, respectively (Table 2 2 ). In site B, the CE and OE were 50% and 25%, respectively, i.e., the user’s and producer’s accuracies were 50% and 75%, respectively (Table 2 2 ). According to the accuracy assessment in these intensive sites, we found that the user’s and producer’s accuracies both ranged from 50% to 75%. Fig. 3 3 shows the partial enlargement of the MODIS-derived rice map and original TM color composite image in which the rice area was obviously labeled out. Visual inspection in Fig. 3 3 revealed significant differences between the rice binary maps derived from MODIS and original TM imagery. Because of the huge differences in spatial resolution and radiometric quality between MODIS and TM data, it is difficult to co-register the TM image with respect to MODIS data by identifying exact common ground control points. Actually, the neighboring pixels around the identified rice distribution region were fatal, resulting in great error when conducting accuracy assessment per pixel. 2 Test field 4 Pixel level 4 Block level (3×3 pixel window) Commission error (%) User’s accuracy (%) Omission error (%) Producer’s accuracy (%) Commission error (%) User’s accuracy (%) Omission error (%) Producer’s accuracy (%) A 26 74 49 51 2 98 21 79 B 50 50 25 75 16 84 6 94 (a, c) MODIS-derived rice maps of sites A and B, respectively; (b, d) Respective enlarged images in the blue rectangle in (a) and (c), respectively; (e, f) Respective TM composite images overlaid by boundary of Rice MOD in (b) and (d). Rice MOD and Rice TM are the rice pixels derived from MODIS and TM imageries, respectively Therefore, in order to compensate for the error caused by geometrical mismatch, we calculated the CE and OE for each 3×3 pixel moving window (Zhan et al., 2000 2000 ). If a pixel was identified as paddy rice with the MODIS algorithm and any pixel in the 3×3 pixel window around it was not labeled as rice in the rice reference map, we consider it as a committed pixel. Otherwise, we consider that the identified pixel based on our algorithm was correct if there was any pixel labeled as rice in the rice reference map. For the calculation of OE, if a pixel was labeled as rice in the rice reference map and any pixel in the 3×3 window around it was not identified as a rice pixel with the MODIS algorithm, we considered it as an omitted pixel. With this criterion, the CE and OE at the 3×3 pixel moving window in site A were 2% and 21%, respectively, while the user’s and producer’s accuracies were 98% and 79%, respectively. In site B, the CE and OE at the 3×3 pixel moving window were 16% and 6%, respectively, while the user’s and producer’s accuracies were 84% and 94%, respectively. Compared with the accuracy at pixel level, the errors were reduced and were more reasonable for evaluating the performance of our algorithm by depressing the spatial mismatching phenomena. The results demonstrated that the MODIS algorithm was very consistent with the actual distribution and reliable in extracting the paddy rice area on a regional scale. Tempo-spatial distribution of paddy rice in Northeast China 2 The annual spatial distribution of paddy rice in Northeast China in the period from 2001 to 2009 was mapped using our algorithm (Fig. 4 4 ). Paddy rice was cultivated in the vast majority of Northeast China except for Daxing’anling City at high latitude. Paddy rice fields in Northeast China were mainly concentrated in three major alluvial plains, viz. Sanjiang Plain, Songnen Plain, and Lower Liaohe Plain, and regions along the river basin with abundant irrigated water supply and favorable meteorological conditions for rice cultivation. The Sanjiang Plain is a low alluvial plain located in northeastern Heilongjiang including Jiamusi, Shuangyashan, Jixi, Hegang, and Qitaihe Cities and developed by three major rivers, namely Heilong River, Songhua River, and Ussuri River with an elevation of 40‒90 m. Songnen Plain is formed by the Songhua River and Nen River and lies in the central region of Northeast China covering the majority of Qiqihar, Suihua, Daqing Cities in the southwest of Heilongjiang Province and Harbin City in Jilin Province. The Lower Liaohe Plain is the alluvial plain of the Liao River and is located in the central part of Liaoning Province. The planting area of paddy rice according to agricultural statistical data in Northeast China has been increasing in general each year, except for 2003. The interannual variation trend in MODIS-derived paddy rice in Northeast China was consistent with the agricultural census data (Fig. 5 5 ). The overall census and MODIS-derived paddy rice acreages in Northeast China increased by about 14 610 km 2 and 14 888 km 2 from 2001 to 2009, respectively. The increase in the rice planting area was mainly concentrated in the Sanjiang Plain with the enhancement of the agricultural irrigated infrastructure and mechanization and breeding of salt tolerant rice varieties. Compared with 2001, the growth rate of the total rice acreage at the end of 2009 was 60%, regardless of the agricultural census or MODIS-derived area. It demonstrated that the paddy rice area derived from MODIS data was confirmed as valid for detecting the rice plant area and evaluating the interannual variation on a large scale. There was a significant and sudden decline in paddy rice in northeast of the Sanjiang Plain in Heilongjiang Province and valley along the East Liao River, Yinma River, and Liu River in the center of Jilin Province in 2003 (Fig. 4 4 ). With the rapid development of rice production, the phenomena of a huge backlog and oversupply of rice had become more serious for farmers in Northeast China in 2002. Furthermore, the decline in the producer price index (PPI) of rice in Liaoning, Jilin and Heilongjiang Provinces in 2002 was a direct economic factor requiring farmers to adjust and optimize the planting structure of paddy rice (Fig. 6 6 ). Data source: China Rural Statistical Yearbook Database ( http://tongji.cnki.net/ http://tongji.cnki.net/ ) Validation of MODIS-derived rice planting area with census data 2 It is a challenge to conduct a field investigation to find the actual planting area of paddy rice in the large extent of Northeast China. Therefore, the published agricultural census data were used to compare the results derived from MODIS data. Fig. 7 7 showed that the total planting area of paddy rice in Northeast China identified using the MODIS-based algorithm has a significant correlation ( R 2 =0.847) with the sown area of paddy rice described in the statistical yearbook from 2001 to 2009. It also indicated that the regression line approached the 1:1 diagonal well. For further analysis, we aggregated and validated the paddy rice area in each city of Liaoning and Jilin Provinces. It should be mentioned that the state farms scattered all over Heilongjiang Province were regarded as an integrated unit in annual provincial agricultural statistics and were not taken into account in the municipal statistical data, despite their geographical location. Thus, the census data at municipal level of Heilongjiang Province were excluded in the further analysis due to its inaccuracy. Fig. 8 8 denoted that the annual paddy rice area derived from the MODIS-based algorithm was in general agreement with municipal statistical data and the scatters concentrated in the vicinity of the 1:1 line. It demonstrated that the high errors in the MODIS-derived area occurred in the cities where the total statistical paddy rice area was less than 300 km 2 and accounted for less than 3% of the territories, such as Benxi, Chaoyang, Dalian, Fushun, and Fuxin Cities in Liaoning Province and Liaoyuan City in Jilin Province. In addition, the rice plant area was overestimated especially in Panjin City, which lies in the central region of the Liaohe Delta and has a high fractional paddy rice plant area, mainly because of the large area of everglades to be found there. Discussion and conclusions 1 An algorithm for extracting a paddy rice field was proposed in this study. The short-wave NIR band (2 105‒2 155 nm) of MOD09A1 was taken into account for the identification of the flooding and transplanting stages of paddy rice with the special spectral characteristics of high water content and low vegetation cover at canopy and satellite level. We compared temporal profiles of MODIS-derived spectral indices of different land cover types at 500-m resolution and hypothesized that a pixel was identified as a potential flooding or transplanting one using the threshold LSWI 2130 > 2EVI. With the validation of paddy rice distribution derived from fine resolution imagery (Landsat TM) at two intensive sites, the user’s and producer’s accuracies were both approximately greater than 80% at the block level (3×3 pixel moving window). The results derived from MODIS data using our improved algorithm were well consistent with the agricultural statistical data on the regional and municipal scales in the period from 2001 to 2009 both in spatial distribution and quantity. This approach is an attempt to highlight the potential of spectral index LSWI 2130 in rice extraction using MODIS time-series analysis. Our algorithm is suitable for analyzing the interannul variability of a rice planting area on the regional scale. Wetland, which is distributed widely in Northeast China, is always difficult to separate from rice. Compared with the results of Sun et al. ( 2009 2009 ), we took the wetland into consideration and improved the ability to exclude wetland. However, there were still deficiencies in the MODIS-derived results with respect to the following factors. Firstly, an MODIS pixel in the MOD09A1 product has an actual spatial resolution of about 463 m and an area of 214 369 m 2 . It is difficult and unreasonable to identify the pixel with a low abundance of paddy rice as a rice pixel. With the exception of the major plains in Northeast China, paddy rice in our study area, especially Jilin Province, often occurs in parcels that are much smaller than a 500-m MODIS pixel, of which the spectral characteristic was an integration of various land cover types. The coarse resolution of MODIS would result in an underestimation of the extent of paddy rice in the regions with complicated topography and fragmentary fields with individual management. Secondly, pixels contaminated by clouds and their shadows in the flooding and transplanting periods over several 8-d periods would be omitted in our algorithm although we had reconstructed the EVI temporal profile to interpolate the contaminated ones. MODIS 8-d composites were generated by selecting the date within the 8-d period with minimum value of the blue band which denoted the clearest atmospheric condition for each pixel. This method might omit some observations associated with the transplanting period. Although the daily reflectance product could be applied to further analysis with its ability for us to obtain the critical flooding and transplanting dates, the trouble with data redundancy and contaminated phenomena is also present and must to be taken into account. Finally, the results might be overestimated in certain regions since a few areas of natural wetland vegetation or other cropland with persistent rainfall in the initial phase might be misclassified as paddy rice, despite typical wetland vegetation having been taken into consideration in this algorithm. Although there are still some uncertainties in the data acquisition and processing, the advantage of large coverage, frequent revisits, and low cost of MODIS data is a better choice for constructing a reliable rice database over the past decade. Farmers always have different rice transplanting schedules, so an algorithm with analysis of time series spectral indices is more suitable for application on a large spatial scale. Our work can also provide a spatial database for estimating yield, net primary production, greenhouse gas emission, and ground water management in a rice field on a large scale. References T Amano , Q Zhu , Y Wang , N Inoue , H Tanaka . Case studies on high yields of paddy rice in Jiangsu Province, China, 1: characteristics of grain production . Jpn J Crop Sci , 1993 . 62 ( 2 ): 267 - 274 . DOI: 10.1626/jcs.62.267 http://doi.org/10.1626/jcs.62.267 . Case studies on high yields of paddy rice in Jiangsu Province, China, 1: characteristics of grain production D Bachelet . Rice paddy inventory in a few provinces of China using AVHRR data . Geocarto Int , 1995 . 10 ( 1 ): 23 - 38 . DOI: 10.1080/10106049509354476 http://doi.org/10.1080/10106049509354476 . Rice paddy inventory in a few provinces of China using AVHRR data JL Barker , JW Burelhach . MODIS Image Simulation from Landsat TM Imagery . ASPRS/ACSM/RT 92 , 1992 . : Fethesda, Washington, DC American Society for Photogrammetry and Remote Sensing , 156 - 165 . . MODIS Image Simulation from Landsat TM Imagery P Ceccato , N Gobron , S Flasse , B Pinty , S Tarantola . Designing a spectral index to estimate vegetation water content from remote sensing data: Part 1. Theoretical approach . Remote Sens Environ , 2002 . 82 ( 2-3 ): 188 - 197 . DOI: 10.1016/S0034-4257(02)00037-8 http://doi.org/10.1016/S0034-4257(02)00037-8 . Designing a spectral index to estimate vegetation water content from remote sensing data: Part 1. Theoretical approach P Ceccato , S Flasse , JM Grégoire . Designing a spectral index to estimate vegetation water content from remote sensing data: Part 2. Validation and applications . Remote Sens Environ , 2002 . 82 ( 2-3 ): 198 - 207 . DOI: 10.1016/S0034-4257(02)00036-6 http://doi.org/10.1016/S0034-4257(02)00036-6 . Designing a spectral index to estimate vegetation water content from remote sensing data: Part 2. Validation and applications K Chandrasekar , MVR Sesha Sai , PS Roy , RS Dwevedi . Land Surface Water Index (LSWI) response to rainfall and NDVI using the MODIS vegetation index product . Int J Remote Sens , 2010 . 31 ( 15 ): 3987 - 4005 . DOI: 10.1080/01431160802575653 http://doi.org/10.1080/01431160802575653 . Land Surface Water Index (LSWI) response to rainfall and NDVI using the MODIS vegetation index product DY Chen , JF Huang , TJ Jackson . Vegetation water content estimation for corn and soybeans using spectral indices derived from MODIS near- and short-wave infrared bands . Remote Sens Environ , 2005 . 98 ( 2-3 ): 225 - 236 . DOI: 10.1016/j.rse.2005.07.008 http://doi.org/10.1016/j.rse.2005.07.008 . Vegetation water content estimation for corn and soybeans using spectral indices derived from MODIS near- and short-wave infrared bands PC Doraiswamy , TR Sinclair , S Hollinger , B Akhmedov , A Stern , J Prueger . Application of MODIS derived parameters for regional crop yield assessment . Remote Sens Environ , 2005 . 97 ( 2 ): 192 - 202 . DOI: 10.1016/j.rse.2005.03.015 http://doi.org/10.1016/j.rse.2005.03.015 . Application of MODIS derived parameters for regional crop yield assessment HL Fang , SL Liang . A hybrid inversion method for mapping leaf area index from MODIS data: experiments and application to broadleaf and needleleaf canopies . Remote Sens Environ , 2005 . 94 ( 3 ): 405 - 424 . DOI: 10.1016/j.rse.2004.11.001 http://doi.org/10.1016/j.rse.2004.11.001 . A hybrid inversion method for mapping leaf area index from MODIS data: experiments and application to broadleaf and needleleaf canopies BC Gao . NDWI–a normalized difference water index for remote sensing of vegetation liquid water from space . Remote Sens Environ , 1996 . 58 ( 3 ): 257 - 266 . DOI: 10.1016/S0034-4257(96)00067-3 http://doi.org/10.1016/S0034-4257(96)00067-3 . NDWI–a normalized difference water index for remote sensing of vegetation liquid water from space MK Gumma , A Nelson , PS Thenkabail , AN Singh . Mapping rice areas of south Asia using MODIS multitemporal data . J Appl Remote Sens , 2011 . 5 ( 1 ): 53547 DOI: 10.1117/1.3619838 http://doi.org/10.1117/1.3619838 . Mapping rice areas of south Asia using MODIS multitemporal data DK Hall , GA Riggs , VV Salomonson , NE DiGirolamo , KJ Bayr . MODIS snow-cover products . Remote Sens Environ , 2002 . 83 ( 1-2 ): 181 - 194 . DOI: 10.1016/S0034-4257(02)00095-0 http://doi.org/10.1016/S0034-4257(02)00095-0 . MODIS snow-cover products AR Huete , HQ Liu , K Batchily , W van Leeuwen . A comparison of vegetation indices over a global set of TM images for EOS-MODIS . Remote Sens Environ , 1997 . 59 ( 3 ): 440 - 451 . DOI: 10.1016/S0034-4257(96)00112-5 http://doi.org/10.1016/S0034-4257(96)00112-5 . A comparison of vegetation indices over a global set of TM images for EOS-MODIS A Huete , K Didan , T Miura , EP Rodriguez , X Gao , LG Ferreira . Overview of the radiometric and biophysical performance of the MODIS vegetation indices . Remote Sens Environ , 2002 . 83 ( 1-2 ): 195 - 213 . DOI: 10.1016/S0034-4257(02)00096-2 http://doi.org/10.1016/S0034-4257(02)00096-2 . Overview of the radiometric and biophysical performance of the MODIS vegetation indices GS Khush . What it will take to feed 5.0 billion rice consumers in 2030 . Plant Mol Biol , 2005 . 59 ( 1 ): 1 - 6 . 16217597 DOI: 10.1007/s11103-005-2159-5 http://doi.org/10.1007/s11103-005-2159-5 . What it will take to feed 5.0 billion rice consumers in 2030 H Kim , AR Huete , PN Nagler , 等 . Monitoring Riparian and Semi-arid Upland Vegetation Using Vegetation and Water Indices from the MODIS Satellite Sensor. Research Insights in Semiarid Ecosystems (RISE) , 2004 . : Tucson University of Arizona , . Monitoring Riparian and Semi-arid Upland Vegetation Using Vegetation and Water Indices from the MODIS Satellite Sensor. Research Insights in Semiarid Ecosystems (RISE) XR Liu , J Zheng , XH Li . Rice purchase and sale market of northeast situation analyze in 2010/2011 . China Grain Econ , 2010 . ( 12 ): 36 - 39 . (in Chinese) . Rice purchase and sale market of northeast situation analyze in 2010/2011 YS Liu , DW Wang , J Gao , W Deng . Land use/cover changes, the environment and water resources in Northeast China . Environ Manage , 2005 . 36 ( 5 ): 691 - 701 . 16206021 DOI: 10.1007/s00267-004-0285-5 http://doi.org/10.1007/s00267-004-0285-5 . Land use/cover changes, the environment and water resources in Northeast China RS Lunetta , JF Knight , J Ediriwickrema , JG Lyon , LD Worthy . Land-cover change detection using multi-temporal MODIS NDVI data . Remote Sens Environ , 2006 . 105 ( 2 ): 142 - 154 . DOI: 10.1016/j.rse.2006.06.018 http://doi.org/10.1016/j.rse.2006.06.018 . Land-cover change detection using multi-temporal MODIS NDVI data L Luo , ZM Wang , KS Song , B Zhang , DW Liu , CY Ren , SM Zhang . Research on the correlation between NDVI and climatic factors of different vegetations in the Northeast China . Acta Bot Bor Occid Sin , 2009 . 29 ( 4 ): 800 - 808 . (in Chinese) . Research on the correlation between NDVI and climatic factors of different vegetations in the Northeast China CL Ma . The study on estimation of vegetation NPP of Xiaerxili nature protection area based on RS . Remote Sens Technol Appl , 2008 . 23 ( 3 ): 323 - 327 . (in Chinese) . The study on estimation of vegetation NPP of Xiaerxili nature protection area based on RS RB Macdonald , FG Hall . Global crop forecasting . Science , 1980 . 208 ( 4445 ): 670 - 679 . 17771086 DOI: 10.1126/science.208.4445.670 http://doi.org/10.1126/science.208.4445.670 . Global crop forecasting JP Malingreau . Global vegetation dynamics: satellite-observations over Asia . Int J Remote Sens , 1986 . 7 ( 9 ): 1121 - 1146 . DOI: 10.1080/01431168608948914 http://doi.org/10.1080/01431168608948914 . Global vegetation dynamics: satellite-observations over Asia HU Neue . Methane emission from rice fields . Bioscience , 1993 . 43 ( 7 ): 466 - 474 . DOI: 10.2307/1311906 http://doi.org/10.2307/1311906 . Methane emission from rice fields K Okamoto , M Fukuhara . Estimation of paddy field area using the area ratio of categories in each mixel of Landsat TM . Int J Remote Sens , 1996 . 17 ( 9 ): 1735 - 1749 . DOI: 10.1080/01431169608948736 http://doi.org/10.1080/01431169608948736 . Estimation of paddy field area using the area ratio of categories in each mixel of Landsat TM DL Peng , AR Huete , JF Huang , FM Wang , HS Sun . Detection and estimation of mixed paddy rice cropping patterns with MODIS data . Int J Appl Earth Obs Geoinf , 2011 . 13 ( 1 ): 13 - 23 . DOI: 10.1016/j.jag.2010.06.001 http://doi.org/10.1016/j.jag.2010.06.001 . Detection and estimation of mixed paddy rice cropping patterns with MODIS data NA Quarmby , M Milnes , TL Hindle , N Silleos . The use of multi-temporal NDVI measurements from AVHRR data for crop yield estimation and prediction . Int J Remote Sens , 1993 . 14 ( 2 ): 199 - 210 . DOI: 10.1080/01431169308904332 http://doi.org/10.1080/01431169308904332 . The use of multi-temporal NDVI measurements from AVHRR data for crop yield estimation and prediction C Rosenzweig , KM Strzepek , DC Major , A Iglesias , DN Yates , A McCluskey , D Hillel . Water resources for agriculture in a changing climate: international case studies . Global Environ Chang , 2004 . 14 ( 4 ): 345 - 360 . DOI: 10.1016/j.gloenvcha.2004.09.003 http://doi.org/10.1016/j.gloenvcha.2004.09.003 . Water resources for agriculture in a changing climate: international case studies DP Roy , JS Borak , S Devadiga , RE Wolfe , M Zheng , J Descloitres . The MODIS land product quality assessment approach . Remote Sens Environ , 2002 . 83 ( 1-2 ): 62 - 76 . DOI: 10.1016/S0034-4257(02)00087-1 http://doi.org/10.1016/S0034-4257(02)00087-1 . The MODIS land product quality assessment approach T Sakamoto , M Yokozawa , H Toritani , M Shibayama , N Ishitsuka , H Ohno . A crop phenology detection method using time-series MODIS data . Remote Sens Environ , 2005 . 96 ( 3-4 ): 366 - 374 . DOI: 10.1016/j.rse.2005.03.008 http://doi.org/10.1016/j.rse.2005.03.008 . A crop phenology detection method using time-series MODIS data DK Sari , IH Ismullah , WN Sulasdi , AB Harto . Detecting rice phenology in paddy fields with complex cropping pattern using time series MODIS data . ITB J Sci , 2010 . 42A ( 2 ): 91 - 106 . DOI: 10.5614/itbj.sci.2010.42.2.2 http://doi.org/10.5614/itbj.sci.2010.42.2.2 . Detecting rice phenology in paddy fields with complex cropping pattern using time series MODIS data Y Shao , XT Fan , H Liu , JH Xiao , S Ross , B Brisco , R Brown , G Staples . Rice monitoring and production estimation using multitemporal RADARSAT . Remote Sens Environ , 2001 . 76 ( 3 ): 310 - 325 . DOI: 10.1016/S0034-4257(00)00212-1 http://doi.org/10.1016/S0034-4257(00)00212-1 . Rice monitoring and production estimation using multitemporal RADARSAT HS Sun , JF Huang , AR Huete , DL Peng , F Zhang . Mapping paddy rice with multi-date moderate-resolution imaging spectroradiometer (MODIS) data in China . J Zhejiang Univ-Sci A , 2009 . 10 ( 10 ): 1509 - 1522 . DOI: 10.1631/jzus.A0820536 http://doi.org/10.1631/jzus.A0820536 . Mapping paddy rice with multi-date moderate-resolution imaging spectroradiometer (MODIS) data in China EF Vermote , A Vermeulen . Atmospheric Correction Algorithm: Spectral Reflectances (MOD09), MODIS Algorithm Technical Background Document, Version 4.0 . 1999 . Department of Geography, University of Maryland , . Atmospheric Correction Algorithm: Spectral Reflectances (MOD09), MODIS Algorithm Technical Background Document, Version 4.0 EF Vermote , N El Saleous , CO Justice , YJ Kaufman , JL Privette , L Remer , JC Roger , D Tanré . Atmospheric correction of visible to middle-infrared EOS-MODIS data over land surfaces: background, operational algorithm and validation . J Geophys Res , 1997 . 102 ( D14 ): 17131 - 17141 . DOI: 10.1029/97JD00201 http://doi.org/10.1029/97JD00201 . Atmospheric correction of visible to middle-infrared EOS-MODIS data over land surfaces: background, operational algorithm and validation Y Wada , W Ohira . Reconstructing Cloud Free SPOT/Vegetation Using Harmonic Analysis with Local Maximum Fitting . 2004 . 25th Asian Conference on Remote Sensing . Tailand 1663 - 1667 . . Reconstructing Cloud Free SPOT/Vegetation Using Harmonic Analysis with Local Maximum Fitting XM Xiao , S Boles , S Frolking , W Salas , B Moore III , C Li , L He , R Zhao . Observation of flooding and rice transplanting of paddy rice fields at the site to landscape scales in China using vegetation sensor data . Int J Remote Sens , 2002 . 23 ( 15 ): 3009 - 3022 . DOI: 10.1080/01431160110107734 http://doi.org/10.1080/01431160110107734 . Observation of flooding and rice transplanting of paddy rice fields at the site to landscape scales in China using vegetation sensor data XM Xiao , S Boles , JY Liu , DF Zhuang , S Frolking , CS Li , W Salas , B Moore . Mapping paddy rice agriculture in southern China using multi-temporal MODIS images . Remote Sens Environ , 2005 . 95 ( 4 ): 480 - 492 . DOI: 10.1016/j.rse.2004.12.009 http://doi.org/10.1016/j.rse.2004.12.009 . Mapping paddy rice agriculture in southern China using multi-temporal MODIS images XM Xiao , S Boles , S Frolking , CS Li , JY Babu , W Salas , B Moore III . Mapping paddy rice agriculture in south and southeast Asia using multi-temporal MODIS images . Remote Sens Environ , 2006 . 100 ( 1 ): 95 - 113 . DOI: 10.1016/j.rse.2005.10.004 http://doi.org/10.1016/j.rse.2005.10.004 . Mapping paddy rice agriculture in south and southeast Asia using multi-temporal MODIS images X Zhan , R Defries , JRG Townshend , C Dimiceli , M Hansen , C Huang , R Sohlberg . The 250 m global land cover change product from the moderate resolution imaging spectroradiometer of NASA’s earth observing system . Int J Remote Sens , 2000 . 21 ( 6-7 ): 1433 - 1460 . DOI: 10.1080/014311600210254 http://doi.org/10.1080/014311600210254 . The 250 m global land cover change product from the moderate resolution imaging spectroradiometer of NASA’s earth observing system LZ Zhang , M Li , ZF Wu , YJ Liu . Vegetation cover change and its mechanism in Northeast China based on SPOT/NDVI data . J Arid Land Resourc Environ , 2011 . 25 ( 1 ): 171 - 175 . (in Chinese) . Vegetation cover change and its mechanism in Northeast China based on SPOT/NDVI data XY Zhang , MA Friedl , CB Schaaf , AH Strahler , JCF Hodges , F Gao , BC Reed , A Huete . Monitoring vegetation phenology using MODIS . Remote Sens Environ , 2003 . 84 ( 3 ): 471 - 475 . DOI: 10.1016/S0034-4257(02)00135-9 http://doi.org/10.1016/S0034-4257(02)00135-9 . Monitoring vegetation phenology using MODIS Y Zhang , CZ Wang , JP Wu , JG Qi , WA Salas . Mapping paddy rice with multitemporal ALOS/PALSAR imagery in southeast China . Int J Remote Sens , 2009 . 30 ( 23 ): 6301 - 6315 . DOI: 10.1080/01431160902842391 http://doi.org/10.1080/01431160902842391 . Mapping paddy rice with multitemporal ALOS/PALSAR imagery in southeast China Y Zhang , YY Wang , SL Su , CS Li . Quantifying methane emissions from rice paddies in Northeast China by integrating remote sensing mapping with a biogeochemical model . Biogeosciences , 2011 . 8 ( 5 ): 1225 - 1235 . DOI: 10.5194/bg-8-1225-2011 http://doi.org/10.5194/bg-8-1225-2011 . Quantifying methane emissions from rice paddies in Northeast China by integrating remote sensing mapping with a biogeochemical model Compliance with ethics guidelines Jing-jing SHI, Jing-feng HUANG, and Feng ZHANG declare that they have no conflict of interest. This article does not contain any studies with human or animal subjects performed by any of the authors. 10.1631/jzus.B1200352 Project supported by the National High-Tech R & D Program (863) of China (No. 2012AA12A30703), the Meteorology Industry Special Project of China Meteorological Administration (CMA) (No. GYHY 201306036), and the Ph.D Programs Foundation of the Ministry of Education of China (No. 20100101110035) 2012-12-19 2013-09-06 2013-06-03 2022-01-05 Copyright © Zhejiang University and Springer-Verlag Berlin Heidelberg 2013 2013 Journal of Zhejiang University-SCIENCE B (Biomedicine and Biotechnology) 10 14