Mapping upto-Date Paddy Rice Extent at 10 M Resolution in China through the Integration of Optical and Synthetic Aperture Radar Images

Rice is a staple food in East Asia and Southeast Asia—an area that accounts for more than half of the world’s population, and 11% of its cultivated land. Studies on rice monitoring can provide direct or indirect information on food security, and water source management. Remote sensing has proven to be the most effective method for the large-scale monitoring of croplands, by using temporary and spectral information. The Google Earth Engine (GEE) is a cloud-based platform providing access to high-performance computing resources for processing extremely large geospatial datasets. In this study, by leveraging the computational power of GEE and a large pool of satellite and other geophysical data (e.g., forest and water extent maps, with high accuracy at 30 m), we generated the first up-to-date rice extent map with crop intensity, at 10 m resolution in the three provinces with the highest rice production in China (the Heilongjiang, Hunan and Guangxi provinces). Optical and synthetic aperture radar (SAR) data were monthly and metric composited to ensure a sufficient amount of up-to-date data without cloud interference. To remove the common confounding noise in the pixel-based classification results at medium to high resolution, we integrated the pixel-based classification (using a random forest classifier) result with the object-based segmentation (using a simple linear iterative clustering (SLIC) method). This integration resulted in the rice planted area data that most closely resembled official statistics. The overall accuracy was approximately 90%, which was validated by ground crop field points. The F scores reached 87.78% in the Heilongjiang Province for monocropped rice, 89.97% and 80.00% in the Hunan Province for monoand double-cropped rice, respectively, and 88.24% in the Guangxi Province for double-cropped rice.


Introduction
Food security is always a major challenge due to continuously increasing populations and limited land availability [1].As a major staple food, rice feeds almost 50% of world's population and provides approximately 20% of the daily human caloric supply due to the high intake of cereal-type foods [2].Monitoring rice cultivation is necessary for understanding the status of food security and providing meaningful information for national food policies and decision-makers [3,4].In addition, paddy rice also plays an important role in water use management because it is the most common water consuming crop.The improvement of water use efficiency in rice planting would greatly help global water security [5,6].Paddy rice areas account for more than 12% of the global cropland area and have spread remarkably since 2002 [7].In 2017, the global rice production (milled) and harvest area were forecasted to reach 500.820 million tons and 165.06 million ha, respectively [8,9] (Figure 1).Asia has the largest number of paddy rice fields and produces more than 90% of the global rice supply [10].China has the largest number of rice planting areas, the greatest grain production, and the greatest rice consumption globally.A total of 65% of people in China use rice as their staple food, such that rice provides more than 50% of the caloric supply [11].
Remote Sens. 2018, 10, x FOR PEER REVIEW 2 of 26 water consuming crop.The improvement of water use efficiency in rice planting would greatly help global water security [5,6].Paddy rice areas account for more than 12% of the global cropland area and have spread remarkably since 2002 [7].In 2017, the global rice production (milled) and harvest area were forecasted to reach 500.820 million tons and 165.06 million ha, respectively [8,9] (Figure 1).Asia has the largest number of paddy rice fields and produces more than 90% of the global rice supply [10].China has the largest number of rice planting areas, the greatest grain production, and the greatest rice consumption globally.A total of 65% of people in China use rice as their staple food, such that rice provides more than 50% of the caloric supply [11].Paddy rice agriculture in China plays a pivotal role in both national and global food security [12].Traditional paddy rice acreage monitoring at a large scale, which requires massive ground surveys and statistical analyses, is time and labor intensive.Remote sensing data can help monitor the ground surface of crops at a large scale by providing precise and timely information on the phenological status and development of vegetation [13,14].Several studies have been conducted on the use of remote sensing data when mapping paddy rice, and the data used by these studies are generally divided into two categories: Optical data, and synthetic aperture radar (SAR) data [15][16][17][18].
Optical data are the most common and easily accessible satellite data, and are often used to explore unique spectral characteristics of a crop at a certain stage using vegetation indices or a combination of indices (e.g., the normalized difference vegetation index (NDVI), the enhanced vegetation index (EVI) [19], and the land surface water index (LSWI) [20].Mapping paddy rice based on optical imagery has a long history that is accompanied by the development of optical satellites [21,22].The unique phenological characteristics of paddy rice (e.g., the transplanting stage using flooding signals) compared to other crops can be used to distinguish other types of natural vegetation.This method is often used to map rice by obtaining either individual images or multiple images during the specific growing stage [23][24][25].It is effective and easy to apply, but obtaining enough optical data without cloud interference during a specific stage at a large scale remains a challenge [26].Furthermore, the time series analysis method tracks the seasonal dynamics of each type of crop using every available image during the crop growing season [27].By using unsupervised clustering methods, several algorithms have been developed to map crops based on their phenological information [28][29][30].The reliance of phenology-based methods on time series and Paddy rice agriculture in China plays a pivotal role in both national and global food security [12].Traditional paddy rice acreage monitoring at a large scale, which requires massive ground surveys and statistical analyses, is time and labor intensive.Remote sensing data can help monitor the ground surface of crops at a large scale by providing precise and timely information on the phenological status and development of vegetation [13,14].Several studies have been conducted on the use of remote sensing data when mapping paddy rice, and the data used by these studies are generally divided into two categories: Optical data, and synthetic aperture radar (SAR) data [15][16][17][18].
Optical data are the most common and easily accessible satellite data, and are often used to explore unique spectral characteristics of a crop at a certain stage using vegetation indices or a combination of indices (e.g., the normalized difference vegetation index (NDVI), the enhanced vegetation index (EVI) [19], and the land surface water index (LSWI) [20].Mapping paddy rice based on optical imagery has a long history that is accompanied by the development of optical satellites [21,22].The unique phenological characteristics of paddy rice (e.g., the transplanting stage using flooding signals) compared to other crops can be used to distinguish other types of natural vegetation.This method is often used to map rice by obtaining either individual images or multiple images during the specific growing stage [23][24][25].It is effective and easy to apply, but obtaining enough optical data without cloud interference during a specific stage at a large scale remains a challenge [26].Furthermore, the time series analysis method tracks the seasonal dynamics of each type of crop using every available image during the crop growing season [27].By using unsupervised clustering methods, several algorithms have been developed to map crops based on their phenological information [28][29][30].The reliance of phenology-based methods on time series and images requires the selection of low to medium spatial and high temporal resolution data [31], which is less capable of identifying crops in small-and medium-sized fields.In China, this limitation is more problematic due to the complexity of small heterogeneous and fragmented farmlands [32].
As a water-intensive crop, paddy rice is usually planted in rainy areas where optical data are highly sensitive to clouds [16].SAR satellite signals can penetrate through clouds and, thus, are suitable for monitoring rice in regions dominated by cloudy and rainy weather.The backscatter in rice fields was significantly lower than that from any other agricultural crop due to the specific water flooding stage during the transplanting stage [33,34].Both single and multi-polarized SAR data have been used to map paddy rice [35,36].However, microwave or radar data are usually expensive and difficult to obtain.Several studies have studied the mapping of rice with SAR data but it still difficult to implement at large areas [37][38][39].
Until recently, monitoring crop dynamics at medium to high resolution has been limited by the lack of available high temporal and spatial resolution satellite imagery [40].The combination of different satellite data will provide increased opportunities for more frequent cloud-free surface observations [41].As the most widely accessible medium to high spatial resolution optical satellite products, the Sentinel-2 and Landsat products are concurrently used to create merged products with high temporal resolutions [42,43].SAR data availability has also increased since the launch of the Sentinel-1A and B on 3 April 2014, and 25 April 2016, respectively.This was the first operational SAR mission that operated within the European Commission's Copernicus program; it provided an unprecedentedly large amount of free data for the operational needs of the Copernicus program [44].The satellite was designed for continuous near real-time land monitoring and provided dual-polarized vertical receiving/vertical transmission with horizontal receiving (VV/VH) SAR images with a global spatial resolution of 5 m to 20 m at least every 10 days.
Processing multisource remote sensing data at high spatial and temporal resolutions on a national scale is a major challenge for both computer storage and operating capacity.High-performance computing systems have become abundant, and large-scale cloud computations have become a universally available commodity [45,46].The Google Earth Engine (GEE) is a cloud-based platform for planetary-scale geospatial analysis that uses Google's massive computational capabilities to study a variety of high-impact societal issues, including deforestation, drought, disaster, disease, food security, water management, climate monitoring, and environmental protection [47].Google has collected a large amount of publicly available remote sensing satellite data from around the world, providing researchers with the ability to apply multisource remote sensing data to various research topics at various spatial scales [48,49].
The objective of this study is to design a universal rice mapping method using the supercomputing power of GEE and massive data sets to generate updated and accurate rice paddy maps with 10 m spatial resolution with crop intensity for China.First, we merged two widely used optical satellite datasets, including Landsat 8 and Sentinel-2, and removed the incidence angle impact in the overlapping areas of the SAR data.Second, we trained a random forest (RF) classifier based on two different types of composites: A phenological monthly composite and a metric composite.Third, object-based segmentations from simple linear iterative clustering (SLIC) were introduced to improve the pixel-based classification.Finally, we quantified the accuracy of the paddy rice map with ground observations and existing algorithms to distinguish rice from other crop observations.This result demonstrates how cloud-based service such as GEE are changing the paradigm of crop monitoring from a static, product-based approach by simplifying data access and processing large amounts of satellite data.

Study Area
China has approximately 120 million ha of agricultural land, of which 25-30% are used for paddy rice production [8].Paddy rice cultivation is practiced in almost every province, the main rice producing areas are in the Northeast, Yangtze River Basin and Pearl river Basin (Figure 2).The number of rice harvests per year gradually increases from one time in the northern provinces up to three times in the southern provinces (rice cultivated by tribes can be planted in parts of the south, but in recent years, these farmers have opted to grow other cash crops instead of rice).Based on the different rice cropping seasons, paddy rice in China can be divided into three typical areas: Mono-, double-cropped rice and mixed.The crop season is determined by the agricultural climate conditions of the planting area.In general, north of the Yangtze River and in the western-most places, only one rice harvest occurs per year.South of the Yangtze Basin, there is a transition zone, with one long or two short rice seasons, which is referred to as the mixed area in this paper.Further south, near the Nanling Mountains, there are generally double crop seasons per year due to the abundance of rainfall and high annual average temperature.In this paper, three provinces, which are the largest producers of the three different rice area types, were selected as the study area: Heilongjiang (monocropped area with 46.00% cloud cover), Hunan (mixed area with 62.39% cloud cover) and Guangxi (double-cropped area with 65.43% cloud cover).

Data and Preprocessing
By taking advantage of the tremendous amounts of petabyte-scale storage for archived remote sensing and related production data as well as the large-scale cloud-based computing capacity of GEE, multisource data in large regions can be chosen and processed in an acceptable amount of time [47].In this study, two widely accessible medium to high spatial resolution multispectral optical datasets, the Landsat 8 Operational Land Imager (OLI) and the Sentinel-2 A/B multispectral instrument (MSI), were chosen as the optical data sources.The new free open-access radar dataset of Sentinel-1 C-band was selected as the SAR data source.All processing and pixeled-based classifications were run in GEE.Based on the different rice cropping seasons, paddy rice in China can be divided into three typical areas: Mono-, double-cropped rice and mixed.The crop season is determined by the agricultural climate conditions of the planting area.In general, north of the Yangtze River and in the western-most places, only one rice harvest occurs per year.South of the Yangtze Basin, there is a transition zone, with one long or two short rice seasons, which is referred to as the mixed area in this paper.Further south, near the Nanling Mountains, there are generally double crop seasons per year due to the abundance of rainfall and high annual average temperature.In this paper, three provinces, which are the largest producers of the three different rice area types, were selected as the study area: Heilongjiang (monocropped area with 46.00% cloud cover), Hunan (mixed area with 62.39% cloud cover) and Guangxi (double-cropped area with 65.43% cloud cover).

Data and Preprocessing
By taking advantage of the tremendous amounts of petabyte-scale storage for archived remote sensing and related production data as well as the large-scale cloud-based computing capacity of GEE, multisource data in large regions can be chosen and processed in an acceptable amount of time [47].In this study, two widely accessible medium to high spatial resolution multispectral optical datasets, the Landsat 8 Operational Land Imager (OLI) and the Sentinel-2 A/B multispectral instrument (MSI), were chosen as the optical data sources.The new free open-access radar dataset of Sentinel-1 C-band was selected as the SAR data source.All processing and pixeled-based classifications were run in GEE.

Cloud-Free Optical Imagery Composition
Currently, Landsat from USGS and Sentinel-2 A/B from the ESA (The European Space Agency) are the most widely accessible medium-high spatial resolution optical products.The synergistic use of these two sources significantly increase the opportunities for timely and cloud-free observations of the surface [43].Blue, green, red, near-infrared (NIR), and two short-wave infrared spectral range (SWIR) bands, which are available in both sensors and have similar spectral response functions (SRF), were used in this paper (Figure 3, Table 1).Currently, Landsat from USGS and Sentinel-2 A/B from the ESA (The European Space Agency) are the most widely accessible medium-high spatial resolution optical products.The synergistic use of these two sources significantly increase the opportunities for timely and cloud-free observations of the surface [43].Blue, green, red, near-infrared (NIR), and two short-wave infrared spectral range (SWIR) bands, which are available in both sensors and have similar spectral response functions (SRF), were used in this paper (Figure 3, Table 1).Land Imager (OLI) and the Sentinel-2 A/B multispectral instrument (MSI).Blue, green, red, nearinfrared (NIR), and two short-wave infrared spectral range (SWIR) bands in each of the two sensors were used in this paper.Detailed spectral information is shown in Table 1.
Table 1.Characteristics of the Sentinel-1 C-band, Sentinel-2 MSI, Landsat-8 OLI band used in this study.Six bands in the optical sensor and two bands in the microwave sensor were used in this study.The elevation and slope layer from the Shuttle Radar Topography Mission (SRTM) were also used for feature extraction.The Hansen Global Forest Change and Joint Research Center (JRC) Global Surface Water Mapping data at a 30 m resolution from the Google Earth Engine (GEE) were used for forest and water mapping in this study.At a pixel level, cloud masking based on a combination of brightness, temperature, and the normalized difference snow index (NDSI) were used to remove clouds from Landsat 8 OLI imagery [50].The QA60 band was used to remove clouds from the Sentinel-2 MSI imagery [51].All data were resampled to 30-m using the average resolution value from all involved pixels.As a result, a total of 8255 images in 269 footprints (7107 from Sentinel-2 and 1148 from Landsat 8) covering the entire paddy rice growing season (February to November) in 2017 were queried from the GEE data pool and used in this study (Table 2).The cloud cover percentage distribution for each province is shown in Figure 4. Due to the rainy weather condition in southern China, the Hunan and Guangxi provinces have greater cloud cover than that in Heilongjiang.At a pixel level, cloud masking based on a combination of brightness, temperature, and the normalized difference snow index (NDSI) were used to remove clouds from Landsat 8 OLI imagery [50].The QA60 band was used to remove clouds from the Sentinel-2 MSI imagery [51].All data were resampled to 30-m using the average resolution value from all involved pixels.As a result, a total of 8255 images in 269 footprints (7107 from Sentinel-2 and 1148 from Landsat 8) covering the entire paddy rice growing season (February to November) in 2017 were queried from the GEE data pool and used in this study (Table 2).The cloud cover percentage distribution for each province is shown in Figure 4. Due to the rainy weather condition in southern China, the Hunan and Guangxi provinces have greater cloud cover than that in Heilongjiang.Collectively, three commonly used indices were calculated and investigated as variables for rice identification.The EVI is a common and useful vegetation index designed to enhance the vegetation signal with an improved sensitivity to soil brightness [52].The LSWI and NDWI are known to be sensitive to the total amounts of liquid water in vegetation and has been widely used to identify paddy rice with the EVI [20].The green chlorophyll vegetation index (GCVI) is typically used to forecast crop yield and has been shown to have a greater dynamic range for denser canopies than that of other vegetation indices.The three indices were defined with the following equations:

Synthetic Aperture Radar (SAR) Image Preprocessing
The Level 1 Ground Range Detected (GRD) product from Sentinel-1 A/B in the IW (Interferometric Wide) swath model used in this study, which has a dual-polarized vertical transmission with VV (Single co-polarization, vertical transmit/vertical receive) and VH (Dual-band cross-polarization, vertical transmit/horizontal receive) bands, had a resolution of 10 m and a swath of 250 km.The incidence angle ranged from approximately 30° to 45°.Each tile had high geometric accuracy and was processed to derive the backscatter coefficient (σ°) in decibels (dB) with GEE (as implemented by the Sentinel-1 Toolbox [53]):

•
Orbit file application (using restituted orbits) Terrain correction (orthorectification) using SRTM 30 or ASTER DEM for areas greater than ±60° latitude, where SRTM was not available.Collectively, three commonly used indices were calculated and investigated as variables for rice identification.The EVI is a common and useful vegetation index designed to enhance the vegetation signal with an improved sensitivity to soil brightness [52].The LSWI and NDWI are known to be sensitive to the total amounts of liquid water in vegetation and has been widely used to identify paddy rice with the EVI [20].The green chlorophyll vegetation index (GCVI) is typically used to forecast crop yield and has been shown to have a greater dynamic range for denser canopies than that of other vegetation indices.The three indices were defined with the following equations:

Synthetic Aperture Radar (SAR) Image Preprocessing
The Level 1 Ground Range Detected (GRD) product from Sentinel-1 A/B in the IW (Interferometric Wide) swath model used in this study, which has a dual-polarized vertical transmission with VV (Single co-polarization, vertical transmit/vertical receive) and VH (Dual-band cross-polarization, vertical transmit/horizontal receive) bands, had a resolution of 10 m and a swath of 250 km.The incidence angle ranged from approximately 30 • to 45 • .Each tile had high geometric accuracy and was processed to derive the backscatter coefficient (σ • ) in decibels (dB) with GEE (as implemented by the Sentinel-1 Toolbox [53]):

•
Orbit file application (using restituted orbits) Terrain correction (orthorectification) using SRTM 30 or ASTER DEM for areas greater than ±60 • latitude, where SRTM was not available.
For this study, the Level 1 GRD product acquired data during the rice crop growing cycle from March to November in 2017.The detailed information for each province is shown in Table 2.
A literature review suggests that VH polarization data are sensitive to rice cultivation [17,54].However, the VH/VV ratio can reduce double-bounce effects and environmental factors (e.g., those due to variations in soil moisture), and may be a more stable indicator than VH or VV backscatter alone [40].In this paper, the VH and VH/VV bands were used for the monthly and metric composites, respectively.
Because the data acquired by GEE were acquired at the province scale, overlapping areas inevitably appeared, as shown in Figure 5a in the yellow area.For the optical data, overlapping areas can significantly increase the number of observations in the time series analysis, thereby improving the quality of information.However, for SAR data, overlapping areas with various incidence angles (approximately 10 • ) could produce noise in time series analysis (Figure 5b).Therefore, incidence angle processing in overlapping areas was conducted before the monthly and metric composites, and the data in overlapping areas with higher incidence angles were removed.Then, a Savitzky-Golay (SG) filter was applied on the temporal axis for each pixel to smooth and reduce system errors.The window size was set as 10, the order and polynomial degree were set as 2 in SG filter (Figure 5c).For this study, the Level 1 GRD product acquired data during the rice crop growing cycle from March to November in 2017.The detailed information for each province is shown in Table 2.
A literature review suggests that VH polarization data are sensitive to rice cultivation [17,54].However, the VH/VV ratio can reduce double-bounce effects and environmental factors (e.g., those due to variations in soil moisture), and may be a more stable indicator than VH or VV backscatter alone [40].In this paper, the VH and VH/VV bands were used for the monthly and metric composites, respectively.
Because the data acquired by GEE were acquired at the province scale, overlapping areas inevitably appeared, as shown in Figure 5a in the yellow area.For the optical data, overlapping areas can significantly increase the number of observations in the time series analysis, thereby improving the quality of information.However, for SAR data, overlapping areas with various incidence angles (approximately 10°) could produce noise in time series analysis (Figure 5b).Therefore, incidence angle processing in overlapping areas was conducted before the monthly and metric composites, and the data in overlapping areas with higher incidence angles were removed.Then, a Savitzky-Golay (SG) filter was applied on the temporal axis for each pixel to smooth and reduce system errors.The window size was set as 10, the order and polynomial degree were set as 2 in SG filter (Figure 5c).  2. (b) Incidence angles of a pixel, which vary greatly in the yellow area where the images overlap.(c) To remove the impact of changing angles in the time series analysis, lower value incidence angles were selected.A Savitzky-Golay (SG) filter was used to smooth the data.

Auxiliary Data
Topographic Data Topographic information such as elevation and slope can provide meaningful data for land cover recognition [55].Generally, slopes of less than 3° are suitable for planting rice paddies due to the effects of flooding during the transplanting phase [56].SRTM 1 arc-second global data at a resolution of 30 m from GEE and its derived variable (slope) were used in this study [57].Figure 6 shows a violin plot of the changes in elevation and slope for the five classes in this study, including water, forest, urban, other crops and rice.The slope of the rice sample mostly lies between 0-5°.(c) To remove the impact of changing angles in the time series analysis, lower value incidence angles were selected.A Savitzky-Golay (SG) filter was used to smooth the data.

Auxiliary Data
Topographic Data Topographic information such as elevation and slope can provide meaningful data for land cover recognition [55].Generally, slopes of less than 3 • are suitable for planting rice paddies due to the effects of flooding during the transplanting phase [56].SRTM 1 arc-second global data at a resolution of 30 m from GEE and its derived variable (slope) were used in this study [57].Figure 6 shows a violin plot of the changes in elevation and slope for the five classes in this study, including water, forest, urban, other crops and rice.The slope of the rice sample mostly lies between 0-5 • .

Stratified Random Sample Points
When training an RF classifier, balanced points from each area and type are necessary.To maintain balanced sample training points, stratified random sampling based on location (the latitude and longitude) was used to obtain forest, water, build-up and paddy rice training samples (the green points in Figure 7).In each province, there are 200 points for each layer were randomly selected for the RF classifier training.The following land cover map was used to sample the training point.

•
Global forest change 2000-2017, version 1.5 The results from the time series analysis of Landsat images characterize the global forest extents and changes from 2000 through 2017.Using version 1.5, the 2017 forest extent was updated, and the Landsat 8 OLI data were used as the data source in recent years.The smallholder rotational agricultural clearing when detecting dry and humid tropical forests was improved.The 2017 forest extent map calculated by tree canopy cover in 2000 and the gains and losses in forest cover in 2017 were used as the forest type in this study [48].

• JRC yearly water classification
The datasets that were available at the time of study are intended to show different facets of the spatial and temporal distributions of surface waters over the last 32 years.The permanent water layer was used as the water type in this study [49].

•
Urban and bare areas Urban and bare areas have either an impervious surface or exposed soils, respectively, which usually have low LSWI values.A frequency map of LSWI < 0 was generated using Equations ( 4) and (5).A pixel with a frequency value > 50% was then identified as either an urban or a bare area [58].
where Flood represents the status of the flooding/transplanting phase, and t represents the period when the observations were acquired.The signals can last for approximately two months after the transplanting phase.The flooding signal that appears during these two months was identified as potential paddy rice.In addition, thresholds for the maximum EVI (>0.60) and minimum EVI (<0.4) during the rice growing season were used to mask sparse and natural vegetation [12].A dynamic range backscatter VH threshold was also used to identify potential rice.The VH sensitivity threshold, which was based on dry reference (P05) and wet reference (P95) images above 9 dB, was used to identify potential rice in this study [16].

Field Data
A total of 2502 field points in Heilongjiang, 606 field points in the Hunan Province and 712 field points in the Guangxi Province, including dry farmland(other crops in accuracy assignment) and paddy rice areas, were collected using GVG (Red points in Figure 7), which is a volunteer agricultural information collection application via smartphone [60,61].Half of the field data were used for training

•
Global forest change 2000-2017, version 1.5 The results from the time series analysis of Landsat images characterize the global forest extents and changes from 2000 through 2017.Using version 1.5, the 2017 forest extent was updated, and the Landsat 8 OLI data were used as the data source in recent years.The smallholder rotational agricultural clearing when detecting dry and humid tropical forests was improved.The 2017 forest extent map calculated by tree canopy cover in 2000 and the gains and losses in forest cover in 2017 were used as the forest type in this study [48].

• JRC yearly water classification
The datasets that were available at the time of study are intended to show different facets of the spatial and temporal distributions of surface waters over the last 32 years.The permanent water layer was used as the water type in this study [49].

•
Urban and bare areas Urban and bare areas have either an impervious surface or exposed soils, respectively, which usually have low LSWI values.A frequency map of LSWI < 0 was generated using Equations ( 4) and (5).A pixel with a frequency value > 50% was then identified as either an urban or a bare area [58]. •

Potential rice map with optical data
The phenological and pixel-based paddy rice mapping (PPPM) method identifies flooding signals during the rice transplanting phase and has been effectively applied and proven in tropical and cold regions.The relationship between LSWI and NDVI (EVI) can effectively identify the flooding and transplanting signals [56,59].Here, we recognize flooding signals using the flowing equation: where Flood represents the status of the flooding/transplanting phase, and t represents the period when the observations were acquired.The signals can last for approximately two months after the transplanting phase.The flooding signal that appears during these two months was identified as potential paddy rice.In addition, thresholds for the maximum EVI (>0.60) and minimum EVI (<0.4) during the rice growing season were used to mask sparse and natural vegetation [12].A dynamic range backscatter VH threshold was also used to identify potential rice.The VH sensitivity threshold, which was based on dry reference (P05) and wet reference (P95) images above 9 dB, was used to identify potential rice in this study [16].

Field Data
A total of 2502 field points in Heilongjiang, 606 field points in the Hunan Province and 712 field points in the Guangxi Province, including dry farmland(other crops in accuracy assignment) and paddy rice areas, were collected using GVG (Red points in Figure 7), which is a volunteer agricultural information collection application via smartphone [60,61].Half of the field data were used for training in the other crop and rice categories.The other half of data including 1910 points were used for the accuracy assessment.

Methodology
A comprehensive overview of the methodology is shown in Figure 8. Elevation data, optical data by harmonized Landsat 8 OLI and Sentinel-2 MSI data and SAR data were integrated as features of the classifier input, and monthly and metric composites were used to ensure the entire geometry; temporary and spectral information were trained in the classifier.An RF classifier was used to generate pixel-based classifications, including forest, water, urban, other crops and paddy rice.A stratified sampling method was used on the forest, water and urban layers to ensure sufficient and balanced samples to train the model and to specifically remove noise confusion, which will be mentioned in the following sections regarding the pixel-based classification results for high-and medium-to high-resolution data.We merged the pixel-based classification results and the object-based segmentation results to remove the impact of noise [62,63].Finally, an accuracy assessment was made based on ground observations and statistical data.stratified sampling method was used on the forest, water and urban layers to ensure sufficient and balanced samples to train the model and to specifically remove noise confusion, which will be mentioned in the following sections regarding the pixel-based classification results for high-and medium-to high-resolution data.We merged the pixel-based classification results and the objectbased segmentation results to remove the impact of noise [62,63].Finally, an accuracy assessment was made based on ground observations and statistical data.

Monthly and Metric Composites
In most regions of China, rice has similar phenological and spectral information to that of other crops, especially maize [32].Normally, only finer temporal resolution time series data can be used to distinguish between rice and other crops.However, due to the cloudy and rainy weather in rice dominate areas, it difficult to collect enough cloud-free optical images for monthly composite at 10 to 30-meter resolution.On the other hand, the phenological information of rice planting in China

Monthly and Metric Composites
In most regions of China, rice has similar phenological and spectral information to that of other crops, especially maize [32].Normally, only finer temporal resolution time series data can be used to distinguish between rice and other crops.However, due to the cloudy and rainy weather in rice dominate areas, it difficult to collect enough cloud-free optical images for monthly composite at 10 to 30-meter resolution.On the other hand, the phenological information of rice planting in China changes every year.Therefore, a metric composite method was used which was invented to capture crop and nature vegetation phenological without the explicit hypothesis or known phenology information of the timing of such dynamics.This composites method is suitable in large areas that it doesn't need to change parameters with location.In this paper, the monthly composites and metric composites were used to extract the object features based on GEE (Figure 9).The strategy was aimed for composition of SAR data and metric composition was based on optical and SAR data.crop and nature vegetation phenological without the explicit hypothesis or known phenology information of the timing of such dynamics.This composites method is suitable in large areas that it doesn't need to change parameters with location.In this paper, the monthly composites and metric composites were used to extract the object features based on GEE (Figure 9).The strategy was aimed for composition of SAR data and metric composition was based on optical and SAR data.Monthly Composites: To detect the time serious phenological information of the paddy rice, monthly composites were used on that.We grouped all processed SAR images into monthly groups from March to November to cover the entire rice growing season in China (9 groups).By leveraging the GEE array-based computational approach, we reduced each monthly group of images into a single median value on a per-pixel, per band basis.The composites results have VH and VH/VV bands with 10 m resolution.
Metric Composites: We generated a series of quantile-based composites for all the bands of optical and SAR imagery by using observations across the entire temporal range.5%, 25%, 50%, 75%, and 95% quantile composites was produced, including the LSWI, EVI and GCVI, leading to a total of 55 metrics (six optical bands, two SAR bands and three indicators).We converted the entire collection of images into a multidimensional array by leveraging the built-in array methods of GEE.By using a percentile reducer, the entire array of observations along the temporal axis (on a per-pixel, per band basis) were sorted, and the quantile values were extracted.Finally, we re-projected the observations corresponding to each quantile back into the images (one per band and quantile).

Pixel-Based Classifier: Random Forest (RF)
The RF classifier uses tree bagging to form an ensemble of trees by searching random subspaces in the given features, then splitting the nodes by minimizing the correlation between the trees [64,65].This method is more robust, provides a relatively faster classification, and is easier to implement than many other classifiers [66].Accurate land cover classifications and improved performances of the RF Monthly Composites: To detect the time serious phenological information of the paddy rice, monthly composites were used on that.We grouped all processed SAR images into monthly groups from March to November to cover the entire rice growing season in China (9 groups).By leveraging the GEE array-based computational approach, we reduced each monthly group of images into a single median value on a per-pixel, per band basis.The composites results have VH and VH/VV bands with 10 m resolution.
Metric Composites: We generated a series of quantile-based composites for all the bands of optical and SAR imagery by using observations across the entire temporal range.5%, 25%, 50%, 75%, and 95% quantile composites was produced, including the LSWI, EVI and GCVI, leading to a total of 55 metrics (six optical bands, two SAR bands and three indicators).We converted the entire collection of images into a multidimensional array by leveraging the built-in array methods of GEE.By using a percentile reducer, the entire array of observations along the temporal axis (on a per-pixel, per band basis) were sorted, and the quantile values were extracted.Finally, we re-projected the observations corresponding to each quantile back into the images (one per band and quantile).

Pixel-Based Classifier: Random Forest (RF)
The RF classifier uses tree bagging to form an ensemble of trees by searching random subspaces in the given features, then splitting the nodes by minimizing the correlation between the trees [64,65].This method is more robust, provides a relatively faster classification, and is easier to implement than many other classifiers [66].Accurate land cover classifications and improved performances of the RF models have been described by many researchers [67][68][69].RF was used to classify rice in this study by using 300 trees with a minimum leaf sample size of 5 to improve classification results.The pixel-based RF classifiers were run with GEE.We trained a classifier for each province, for each category there are 200 points comes from the stratified random samples and for the other crops and rice, there are also fields points have been added.

Simple Linear Iterative Clustering (SLIC) Superpixel Segmentation
Image segmentation is a principle function that splits an image into separated regions or objects depending on the specified parameters.Pixels with similar spectral and spatial value are considered an object according to the object-based classifier.The segmentation technique utilizes spatial concepts that involve geometric features, spatial relations, and scaled topological relations between upscale and downscale inheritances [70,71].Classical segmentation algorithms mainly include the recursive hierarchical segmentation (RHSeg) algorithm [72,73], multiresolution segmentation [74] and watershed segmentation [75].The SLIC clusters pixels in combined five-dimensional color and the image plane space to efficiently generate compact and nearly uniform superpixels.This method is fast and computationally efficient because the SLIC does not compare each pixel with all pixels in the scene [76].In this paper, an SLIC segmentation was used to develop field segmentations [77].The true color red, green and blue bands from the Sentinel-2 MSI imagery at a 10 m resolution were used and stretched to 0-255.The compactness was set to 20, and sigma was set to 1.5 to balance color and spatial proximity and to smooth the edge of each segmented area.The raw data and the segmented field results are shown in Figure 10a,c.

Integration of the Pixel-Based Classification and the Object-Based Segmentation
The pixel-based clustering algorithms focused only on the spectral value of each pixel, which resulted in different types of noise confusion when applied to a high-resolution image; this occurrence has been referred to as "salt and pepper" noise [78].In this paper, we integrated the pixel-based classification and the object-based segmentation to improve the classification and remove the "salt and pepper" noise.Each segmented area in the output of the SLIC consisted of pixels with a unique label, which were labeled "rice" or "no rice" when the segment was covered with the pixel-based classification results.To merge the pixel-based classification and the field boundary from the segmentation, the pixel-based classification results were set to 0 (no rice) and 1 (rice).The mean value of each segmented area was calculated by the pixel-based classification, and a threshold was set to ensure that each segmented area was labeled [79]: label = mean ≥ 0.6, rice mean < 0.6, no rice (7) The example shown in Figure 10 highlights the resulting value after the merging process.The center coordinates of the example area (the Heilongjiang Province) are 132 • 49 E and 47 • 03 N, where square-shaped fields usually indicate rice, and fields in long strips usually indicate maize.Maize is primarily planted on both sides of the road, as shown in the central part of the image (Figure 10a).The pixel-based classification of rice (blue) covered most rice areas.However, some pixels were missing, and some appeared in the maize area due to cropland heterogeneity and spectral contamination among neighboring pixels (Figure 10b).Figure 10d was merged to produce more refined and complete boundaries of the rice fields (blue), which improved in consistency.
image plane space to efficiently generate compact and nearly uniform superpixels.This method is fast and computationally efficient because the SLIC does not compare each pixel with all pixels in the scene [76].In this paper, an SLIC segmentation was used to develop field segmentations [77].The true color red, green and blue bands from the Sentinel-2 MSI imagery at a 10 m resolution were used and stretched to 0-255.The compactness was set to 20, and sigma was set to 1.5 to balance color and spatial proximity and to smooth the edge of each segmented area.The raw data and the segmented field results are shown in Figure 10a,c.

Accuracy Assessment
The field points were used for accuracy assignment, there are 1910 points which were defined as the other crop and rice categories, used for the accuracy assessment.Four different metrics, which were derived from the confusion matrix, were selected for accuracy assessment.The overall accuracy (OA) evaluated the overall effectiveness of the algorithm, and the F score measured the accuracy of a class using precision and recall measures.The study established error matrices in each province, which calculated the OA, user accuracies (UA), producer accuracies (PA), and F score using the following equations: where S d represents the total number of correctly classified pixels, n represents the total number of validation pixels, and X ij represents an observation in row i and column j in the confusion matrix; X i represents the marginal total of row i, and X j represents the marginal total of column j in the confusion matrix.

Results
The study produced a nominal 10 m rice extent product with crop intensities in three typical provinces in China by integrating optical and SAR data for the year 2017 (Figure 11).The areas analyzed in this study were compared before and after integration of the pixel-based classification and object-based segmentation, where the area data were reported by the National Bureau of Statistics of China.The rice area based on the pixel-based classification agreed with the statistical rice area.However, the area results were higher than the statistics results.By merging with the segmented field, the misclassed pixels were removed, the result is closer to the statistics value (Figure 12).In the following subsections, we discuss the accuracy of our results in the Heilongjiang, Hunan and Guangxi provinces individually.
Remote Sens. 2018, 10, x FOR PEER REVIEW 15 of 26 area.However, the area results were higher than the statistics results.By merging with the segmented field, the misclassed pixels were removed, the result is closer to the statistics value (Figure 12).In the following subsections, we discuss the accuracy of our results in the Heilongjiang, Hunan and Guangxi provinces individually.

Monthly and Metric Composites for Each Province
Figure 13a,c,e show the differences in the quantiles for each class averaged across all the training points.Compared with the true color bands (i.e., red, green, and blue), NIR, SWIR1, and SWIR2 showed improved diversities for each category due to their higher sensitivities to vegetation.The LSWI showed large differences between rice and other crops.The EVI and GCVI showed adequate diversities between vegetation and other categories.Similarly, the differences in phenology for each category are shown in Figure 13b,d,f, which show VH/VV and VH band time series averaged during all paddy rice growing stages across all the training points.During the transplanting stage of paddy rice (May and June) monocropped rice (April), and double-cropped rice (August), water flooding significantly decreased backscatter, which helped to identify paddy rice and its cropping intensity.

Monthly and Metric Composites for Each Province
Figure 13a,c,e show the differences in the quantiles for each class averaged across all the training points.Compared with the true color bands (i.e., red, green, and blue), NIR, SWIR1, and SWIR2 showed improved diversities for each category due to their higher sensitivities to vegetation.The LSWI showed large differences between rice and other crops.The EVI and GCVI showed adequate diversities between vegetation and other categories.Similarly, the differences in phenology for each category are shown in Figure 13b,d,f

Monthly and Metric Composites for Each Province
Figure 13a,c,e show the differences in the quantiles for each class averaged across all the training points.Compared with the true color bands (i.e., red, green, and blue), NIR, SWIR1, and SWIR2 showed improved diversities for each category due to their higher sensitivities to vegetation.The LSWI showed large differences between rice and other crops.The EVI and GCVI showed adequate diversities between vegetation and other categories.Similarly, the differences in phenology for each category are shown in Figure 13b,d,f, which show VH/VV and VH band time series averaged during all paddy rice growing stages across all the training points.During the transplanting stage of paddy rice (May and June) monocropped rice (April), and double-cropped rice (August), water flooding significantly decreased backscatter, which helped to identify paddy rice and its cropping intensity.

Accuracy Assessment
This rice extent product in the Heilongjiang, Hunan and Guangxi provinces was systematically tested for accuracy (Table 3) based on field data with other crops and rice with crop intensity.There are 1910 points (1251 in Heilongjiang, 303 in Hunan, and 356 in Guangxi) which were defined as the other crop and rice categories, used for the accuracy assessment.
As shown in Figure 13a, rice has similar features in its metric composites.In the Heilongjiang and Guangxi provinces, when operational monocropped or double-cropped rice was planted, the OA were 93.37% and 95.51%, with PA of 90.03% and 92.31%, and UA of 85.63% and 84.51% for rice (Table 3), respectively.In the Hunan Province, monocropped and double-cropped rice were mixed when planted.Generally, monocropped rice is planted north of the Yangtze River, and double-cropped rice is planted south of the river due to climatic conditions.Our results (Figure 11) clearly show this pattern, with OA reaching 89.11%; PA were 86.67% and 89.29% and UA were 93.53% and 72.46% for monocropped and double-cropped rice, respectively.The F scores varied between 89.97% and 80.00%.

Accuracy Assessment
This rice extent product in the Heilongjiang, Hunan and Guangxi provinces was systematically tested for accuracy (Table 3) based on field data with other crops and rice with crop intensity.There are 1910 points (1251 in Heilongjiang, 303 in Hunan, and 356 in Guangxi) which were defined as the other crop and rice categories, used for the accuracy assessment.
As shown in Figure 13a, rice has similar features in its metric composites.In the Heilongjiang and Guangxi provinces, when operational monocropped or double-cropped rice was planted, the OA were 93.37% and 95.51%, with PA of 90.03% and 92.31%, and UA of 85.63% and 84.51% for rice (Table 3), respectively.In the Hunan Province, monocropped and double-cropped rice were mixed when planted.Generally, monocropped rice is planted north of the Yangtze River, and double-cropped rice is planted south of the river due to climatic conditions.Our results (Figure 11) clearly show this pattern, with OA reaching 89.11%; PA were 86.67% and 89.29% and UA were 93.53% and 72.46% for monocropped and double-cropped rice, respectively.The F scores varied between 89.97% and 80.00%.

Discussion
Several studies have globally mapped croplands at a 30 m resolution and determined acceptable accuracies [79,80].However, identifying specific crop types is still challenging.This study indicated the feasibility and reliability of mapping up-to-date annual 10-m paddy rice in three typical and high production provinces by integrating the currently available optical and SAR images.To our knowledge, this study was the first to investigate paddy rice mapping at a 10 m spatial resolution on a national scale (the largest producer provinces of the three different area types were selected as the study area).The feasibility of this study was attributed to several factors, including improved data quality and quantity and used of a cloud-based platform and a machine learning classifier.
Comparisons among the Proba-V-based 100 m resolution rice map (based on the PPPM algorithm, which identifies paddy rice by determining flooding/transplanting signals, has been the most successful and widely used rice monitoring method that utilizes optical images), the pixel-based 30 m rice map with harmonized Landsat 8 and Sentinel-2 images (using the RF classifier) and the merged 10 m rice map (with SLIC segmentation and pixel-based classification) in a typical paddy rice cropping region in Heilongjiang Province (110°14′E and 47°39′N) are shown in Figure 14.Because of 10-30 m high-resolution optical data, the pixel-based 30 m rice map was able to resolve more details with higher accuracy (Figure 14a,b).The pixel-based classification results tended to possess "salt and pepper" noise; that is, several categories exhibited noise characteristics due to spectral confusion among the various land cover classes (Figure 14b with red circle) [78].The merged 10 m rice map removed this noise effectively, especially in urban areas (Figure 14c).
Flooding/transplanting signals, which have been the most important feature used in rice monitoring, can last for approximately two months after transplanting.However, this relative approach based on optical imagery requires a sufficient amount of clear observations during this stage.In southern China, rice is usually planted in areas with rainy weather, which leads to cloudy observations.Figure 15 shows the number of clear observation pixels in Guangxi Province with harmonized Landsat 8 and Sentinel-2 images spanning from the rice transplanting phase to the early vegetative growth phase of March and April.A total of 59.1% of areas did not have clear observation

Discussion
Several studies have globally mapped croplands at a 30 m resolution and determined acceptable accuracies [79,80].However, identifying specific crop types is still challenging.This study indicated the feasibility and reliability of mapping up-to-date annual 10-m paddy rice in three typical and high production provinces by integrating the currently available optical and SAR images.To our knowledge, this study was the first to investigate paddy rice mapping at a 10 m spatial resolution on a national scale (the largest producer provinces of the three different area types were selected as the study area).The feasibility of this study was attributed to several factors, including improved data quality and quantity and used of a cloud-based platform and a machine learning classifier.
Comparisons among the Proba-V-based 100 m resolution rice map (based on the PPPM algorithm, which identifies paddy rice by determining flooding/transplanting signals, has been the most successful and widely used rice monitoring method that utilizes optical images), the pixel-based 30 m rice map with harmonized Landsat 8 and Sentinel-2 images (using the RF classifier) and the merged 10 m rice map (with SLIC segmentation and pixel-based classification) in a typical paddy rice cropping region in Heilongjiang Province (110 • 14 E and 47 • 39 N) are shown in Figure 14.Because of 10-30 m high-resolution optical data, the pixel-based 30 m rice map was able to resolve more details with higher accuracy (Figure 14a,b).The pixel-based classification results tended to possess "salt and pepper" noise; that is, several categories exhibited noise characteristics due to spectral confusion among the various land cover classes (Figure 14b with red circle) [78].The merged 10 m rice map removed this noise effectively, especially in urban areas (Figure 14c).pixels, and only 31.1% of the area had a single observation, making it difficult to identify flooding signatures using optical imagery in these areas.Flooding/transplanting signals, which have been the most important feature used in rice monitoring, can last for approximately two months after transplanting.However, this relative approach based on optical imagery requires a sufficient amount of clear observations during this stage.In southern China, rice is usually planted in areas with rainy weather, which leads to cloudy observations.Figure 15 shows the number of clear observation pixels in Guangxi Province with harmonized Landsat 8 and Sentinel-2 images spanning from the rice transplanting phase to the early vegetative growth phase of March and April.A total of 59.1% of areas did not have clear observation pixels, and only 31.1% of the area had a single observation, making it difficult to identify flooding signatures using optical imagery in these areas.SAR imaging made it possible to map the most current rice data in cloudy areas.Several investigations have demonstrated that the C-band is a useful parameter for mapping and monitoring rice croplands.Nguyen and Wagner set a dynamic range backscatter threshold for images (8.5 dB) to identify potential rice areas in the Mediterranean region [16].A violin plot showing the changes in VH sensitivity (VHsen) based on dry (P05) and wet (P95) references during the rice growing season for 5 classes in the Heilongjiang, Hunan and Guangxi provinces is shown in Figure 16.Due to the paddy rice flooded stage, the VH sensitivity is higher than that of other types as previous studies have shown [17,33,81].The dynamic range backscatter threshold works in small scale benchmarks.However, it remains difficult to apply at the province or bigger scale.The addition of SAR data, especially the monthly time series composite based on backscatter, greatly improves the generalization and adaptability of the method.
GEE cloud-based computing offers substantial computing power by linking all computers in the Google data center, which allows parallel processing and thus enables the classification of provinces with 30 m pixels in a matter of hours.In this study, the majority of the work was accomplished using GEE, in which parts of the study are still running at the local level, including the segmentation and accuracy assessment.The orientation-based classification, based on image segmentation, has been proven to remove the noise that appears in pixel-based classifications with medium to high resolution images.This was also confirmed in this study.Until now, it has been difficult to perform segmentation using GEE directly.
Compared with the pixel-based classification results, the object-oriented classification, which was used by integrating the object-based segmentation, removed the noise between adjacent categories and improved the accuracy in predicting rice areas.However, the satellite imagery with a 10 m resolution and the segmentation method are still limited for small fragmented fields and other objects, especially in southern China, where farmland boundaries rarely exist, and mixed crops are often planted in the same fields.In such cases, SLIC is sensitive to texture and can generate smooth, regular-sized superpixels in non-textured regions and highly irregular superpixels in textured regions [82].We have to adjust the parameters including the desired number to be generated and the compactness of the superpixels in such areas.This method increases human involvement and lacks SAR imaging made it possible to map the most current rice data in cloudy areas.Several investigations have demonstrated that the C-band is a useful parameter for mapping and monitoring rice croplands.Nguyen and Wagner set a dynamic range backscatter threshold for images (8.5 dB) to identify potential rice areas in the Mediterranean region [16].A violin plot showing the changes in VH sensitivity (VHsen) based on dry (P05) and wet (P95) references during the rice growing season for 5 classes in the Heilongjiang, Hunan and Guangxi provinces is shown in Figure 16.Due to the paddy rice flooded stage, the VH sensitivity is higher than that of other types as previous studies have shown [17,33,81].The dynamic range backscatter threshold works in small scale benchmarks.However, it remains difficult to apply at the province or bigger scale.The addition of SAR data, especially the monthly time series composite based on backscatter, greatly improves the generalization and adaptability of the method.
GEE cloud-based computing offers substantial computing power by linking all computers in the Google data center, which allows parallel processing and thus enables the classification of provinces with 30 m pixels in a matter of hours.In this study, the majority of the work was accomplished using GEE, in which parts of the study are still running at the local level, including the segmentation and accuracy assessment.The orientation-based classification, based on image segmentation, has been proven to remove the noise that appears in pixel-based classifications with medium to high resolution images.This was also confirmed in this study.Until now, it has been difficult to perform segmentation using GEE directly.
generalization and adaptability in assessing large areas.Additionally, SLIC segmentation fails to identify the boundaries of crop fields using 10 m imagery.Recently, deep learning has become a technology of considerable interest for computer vision tasks such as object recognition and classification [83] that rely on the rapid development of Graphic Processing Units (GPU).Semantic segmentation-based methods, such as fully convolutional networks (FCNs) [84], SegNet [85], and U-Net [86], have attracted substantial attention due to their efficiency in classification and segmentation of images.These methods effectively and quickly assign a semantic label (i.e., a class) to each coherent region of an image and provide opportunities for remote sensing imagery recognition and mapping.By adapting contemporary classification networks to FCNs and transferring their learned representations by fine-tuning the segmentation task, this Compared with the pixel-based classification results, the object-oriented classification, which was used by integrating the object-based segmentation, removed the noise between adjacent categories and improved the accuracy in predicting rice areas.However, the satellite imagery with a 10 m resolution and the segmentation method are still limited for small fragmented fields and other objects, especially in southern China, where farmland boundaries rarely exist, and mixed crops are often planted in the same fields.In such cases, SLIC is sensitive to texture and can generate smooth, regular-sized superpixels in non-textured regions and highly irregular superpixels in textured regions [82].We have to adjust the parameters including the desired number to be generated and the compactness of the superpixels in such areas.This method increases human involvement and lacks generalization and adaptability in assessing large areas.Additionally, SLIC segmentation fails to identify the boundaries of crop fields using 10 m imagery.
Recently, deep learning has become a technology of considerable interest for computer vision tasks such as object recognition and classification [83] that rely on the rapid development of Graphic Processing Units (GPU).Semantic segmentation-based methods, such as fully convolutional networks (FCNs) [84], SegNet [85], and U-Net [86], have attracted substantial attention due to their efficiency in classification and segmentation of images.These methods effectively and quickly assign a semantic label (i.e., a class) to each coherent region of an image and provide opportunities for remote sensing imagery recognition and mapping.By adapting contemporary classification networks to FCNs and transferring their learned representations by fine-tuning the segmentation task, this method merges segmentation and recognition into one training step.By training with a massive dataset, the model can become much more generalizable and adaptable.A combination of GEE and the TensorFlow deep learning framework is a part of Google's future development [47].We have reason to believe that further work on crop recognition and mapping, not limited to paddy rice, can be performed using the deep learning technology.

Conclusions
The intent of this study was to demonstrate how new cloud-based service, such as GEE, allow users to generate accurate cropping maps of large areas at medium to high resolutions based on multi-source data and machine learning classifiers.By leveraging the computational power of GEE, a large pool of satellites and other geophysical data (e.g., the forest and water extent map at 30 m with a high accuracy), we generated an up-to-date rice extent map of crop intensity at a 10 m resolution in three typical provinces of high production in China for the first time (the Heilongjiang, Hunan and Guangxi provinces).Optical and SAR data were composited monthly and metrically to ensure a sufficient amount of up-to-date data without cloud interference.Based on the integration of pixel-based classifications (using the RF classifier) with the object-based segmentation (using SLIC), commonly observed noise in the pixel-based classification at medium to high resolutions was removed, and the planted rice area more closely resembled that of the official statistics.The OA, which was validated by 1910 ground crop field points, reached approximately 90%.The F scores reached 87.78% in Heilongjiang Province for monocropped rice; 89.97% and 80.00% in Hunan Province for mono-and double-cropped rice, respectively; and 88.24% in Guangxi Province for double-cropped rice.

Figure 1 .
Figure 1.The global rice area and annual change in production from The Food and Agriculture Organization and the Agricultural Market Information Service (FAO-AMIS).

Figure 1 .
Figure 1.The global rice area and annual change in production from The Food and Agriculture Organization and the Agricultural Market Information Service (FAO-AMIS).
Remote Sens. 2018, 10, x FOR PEER REVIEW 4 of 26 number of rice harvests per year gradually increases from one time in the northern provinces up to three times in the southern provinces (rice cultivated by tribes can be planted in parts of the south, but in recent years, these farmers have opted to grow other cash crops instead of rice).

Figure 2 .
Figure 2. The location of the study area (China) and the rice production of each province in 2015 (the production data of each province are from the National Bureau of Statistics in China).Three types of regions were identified based on different cropped seasons: A monocropped region, a doublecropped region and a mixed region.The Heilongjiang, Hunan and Guangxi provinces, which are the largest producers of the three different area types, were selected as the study area in this paper.The right column in the figure shows the annual normalized difference vegetation index (NDVI) growth process line and standard deviation for each rice type.

Figure 2 .
Figure 2. The location of the study area (China) and the rice production of each province in 2015 (the production data of each province are from the National Bureau of Statistics in China).Three types of regions were identified based on different cropped seasons: A monocropped region, a double-cropped region and a mixed region.The Heilongjiang, Hunan and Guangxi provinces, which are the largest producers of the three different area types, were selected as the study area in this paper.The right column in the figure shows the annual normalized difference vegetation index (NDVI) growth process line and standard deviation for each rice type.

Figure 3 .
Figure 3. Average relative spectral response functions (SRF) of bands using the Landsat 8 Operational Land Imager (OLI) and the Sentinel-2 A/B multispectral instrument (MSI).Blue, green, red, nearinfrared (NIR), and two short-wave infrared spectral range (SWIR) bands in each of the two sensors were used in this paper.Detailed spectral information is shown in Table1.

Figure 3 .
Figure 3. Average relative spectral response functions (SRF) of bands using the Landsat 8 Operational Land Imager (OLI) and the Sentinel-2 A/B multispectral instrument (MSI).Blue, green, red, near-infrared (NIR), and two short-wave infrared spectral range (SWIR) bands in each of the two sensors were used in this paper.Detailed spectral information is shown in Table1.

Figure 4 .
Figure 4.The number of observation pixels in the Heilongjiang, Hunan and Guangxi provinces by harmonized Landsat 8 and Sentinel-2 images between March and November of 2017.The detailed and required scene information is shown in Table 2.The right column shows the cloud cover percentage distribution for each image.

Figure 4 .
Figure 4.The number of observation pixels in the Heilongjiang, Hunan and Guangxi provinces by harmonized Landsat 8 and Sentinel-2 images between March and November of 2017.The detailed and required scene information is shown in Table 2.The right column shows the cloud cover percentage distribution for each image.

Figure 5 .
Figure 5. (a) The number of observation pixels using the Sentinel-1 C-band imager in the Hunan Province between March and November of 2017.The detailed required scene information is shown in Table 2. (b) Incidence angles of a pixel, which vary greatly in the yellow area where the images overlap.(c) To remove the impact of changing angles in the time series analysis, lower value incidence angles were selected.A Savitzky-Golay (SG) filter was used to smooth the data.

Figure 5 .
Figure 5. (a) The number of observation pixels using the Sentinel-1 C-band imager in the Hunan Province between March and November of 2017.The detailed required scene information is shown in Table 2. (b) Incidence angles of a pixel, which vary greatly in the yellow area where the images overlap.(c)To remove the impact of changing angles in the time series analysis, lower value incidence angles were selected.A Savitzky-Golay (SG) filter was used to smooth the data.

Figure 6 .
Figure 6.Violin plot of changes in elevation and slope for the five classes analyzed in this study, including water, forest, urban, other crops and rice.Because paddy rice cannot be planted in sloped areas due to the effects of floods during the transplanting phase, the slope of the rice sample was mostly kept between 0-5°.

Figure 6 .•
Figure 6.Violin plot of changes in elevation and slope for the five classes analyzed in this study, including water, forest, urban, other crops and rice.Because paddy rice cannot be planted in sloped areas due to the effects of floods during the transplanting phase, the slope of the rice sample was mostly kept between 0-5 • .

Figure 7 .
Figure 7.The spatial distribution of crop field points and sampling points (green) by stratified sampling from the water, forest and urban distribution maps.The red points represent crop field points that are either dry croplands or paddy rice regions collected using GVG (GPS, Video, and GIS, a volunteer agricultural information collection application via smartphones); 50% of the points were used to train the model, and the other points were used for validation.

Figure 7 .
Figure 7.The spatial distribution of crop field points and sampling points (green) by stratified sampling from the water, forest and urban distribution maps.The red points represent crop field points that are either dry croplands or paddy rice regions collected using GVG (GPS, Video, and GIS, a volunteer agricultural information collection application via smartphones); 50% of the points were used to train the model, and the other points were used for validation.

Figure 8 .
Figure 8. Overview of the methodology for cropland extent mapping.The study integrates a pixelbased classification involving random forests (RFs) with object-oriented simple linear iterative clustering (SLIC).The chart also shows the reference and training dataset used.

Figure 8 .
Figure 8. Overview of the methodology for cropland extent mapping.The study integrates a pixel-based classification involving random forests (RFs) with object-oriented simple linear iterative clustering (SLIC).The chart also shows the reference and training dataset used.

Figure 9 .
Figure 9. Monthly and metric composites by leveraging the GEE array-based computational approach.

Figure 9 .
Figure 9. Monthly and metric composites by leveraging the GEE array-based computational approach.

Figure 10 .
Figure 10.The Heilongjiang (132°49′E and 47°03′N) example of (a) the Sentinel-2 MSI color composite image, which had the lowest cloud cover during the rice growing season of 2017; (b) the pixel-based classification using the random forest (RF) classifier; (c) the object-based SLIC image segmentation result and (d) the merged results, which merged the SLIC segmentation result with the pixel-based RF classification.

Figure 10 .
Figure 10.The Heilongjiang (132 • 49 E and 47 • 03 N) example of (a) the Sentinel-2 MSI color composite image, which had the lowest cloud cover during the rice growing season of 2017; (b) the pixel-based classification using the random forest (RF) classifier; (c) the object-based SLIC image segmentation result and (d) the merged results, which merged the SLIC segmentation result with the pixel-based RF classification.

Figure 11 .
Figure 11.Resultant rice maps of the Heilongjiang, Hunan, and Guangxi provinces with crop intensity.Figures in the right column show local details in different provinces using Sentinel-2 images at a 10 m resolution.

Figure 11 .
Figure 11.Resultant rice maps of the Heilongjiang, Hunan, and Guangxi provinces with crop intensity.Figures in the right column show local details in different provinces using Sentinel-2 images at a 10 m resolution.

Figure 12 .
Figure 12.Scatterplot of the 2015 and 2017 rice harvest area statistics used in this study, which shows that the integration of the pixel-based classification and the object-based segmentation has an improved agreement over statistical data alone (National Bureau of Statistics of China).

Figure 12 .
Figure 12.Scatterplot of the 2015 and 2017 rice harvest area statistics used in this study, which shows that the integration of the pixel-based classification and the object-based segmentation has an improved agreement over statistical data alone (National Bureau of Statistics of China).

26 Figure 12 .
Figure13a,c,e show the differences in the quantiles for each class averaged across all the training points.Compared with the true color bands (i.e., red, green, and blue), NIR, SWIR1, and SWIR2 showed improved diversities for each category due to their higher sensitivities to vegetation.The LSWI showed large differences between rice and other crops.The EVI and GCVI showed adequate diversities between vegetation and other categories.Similarly, the differences in phenology for each category are shown in Figure13b,d,f, which show VH/VV and VH band time series averaged during all paddy rice growing stages across all the training points.During the transplanting stage of paddy rice (May and June) monocropped rice (April), and double-cropped rice (August), water flooding significantly decreased backscatter, which helped to identify paddy rice and its cropping intensity.

Figure 13 .
Figure 13.(a,c,e) shows the quantiles of different spectral properties for each category based on metrics composites, while (b,d,f) shows average monthly SAR spectral information for paddy rice and other crop, as extracted from the monthly composites.Each line shows the mean with their standard deviation from all the training points for each category.

Figure 13 .
Figure 13.(a,c,e) shows the quantiles of different spectral properties for each category based on metrics composites, while (b,d,f) shows average monthly SAR spectral information for paddy rice and other crop, as extracted from the monthly composites.Each line shows the mean with their standard deviation from all the training points for each category.

Figure 14 .
Figure 14.Visual comparisons among the Proba-V-based 100 m rice map (based on the PPPM algorithm, (a)), the pixel-based 30-m rice map (using a random forest (RF) classifier, (b)) and the merged 10 m rice map (using the SLIC segmentation and the pixel-based classification, (c)) for a typical paddy rice cropping region in Heilongjiang Province (110°14′E and 47°39′N).(d) was the Sentinel-2 MSI color composite image with the lowest cloud cover during the 2017 rice growing season.The right column showcases the spatial details of (a-d) within the red box.

Figure 14 .
Figure 14.Visual comparisons among the Proba-V-based 100 m rice map (based on the PPPM algorithm, (a)), the pixel-based 30-m rice map (using a random forest (RF) classifier, (b)) and the merged 10 m rice map (using the SLIC segmentation and the pixel-based classification, (c)) for a typical paddy rice cropping region in Heilongjiang Province (110 • 14 E and 47 • 39 N).(d) was the Sentinel-2 MSI color composite image with the lowest cloud cover during the 2017 rice growing season.The right column showcases the spatial details of (a-d) within the red box.

26 Figure 15 .
Figure 15.The number of clear observation pixels in Guangxi Province with harmonized Landsat 8 and Sentinel-2 images from the rice transplanting phase to the early vegetative growth phase between March and April.A total of 59.1% of areas did not have observation pixels, and only 31.1% of areas had a single observation, making it difficult to find flooding.

Figure 15 .
Figure 15.The number of clear observation pixels in Guangxi Province with harmonized Landsat 8 and Sentinel-2 images from the rice transplanting phase to the early vegetative growth phase between March and April.A total of 59.1% of areas did not have observation pixels, and only 31.1% of areas had a single observation, making it difficult to find flooding.

Figure 16 .
Figure 16.Violin plot of changes in VH sensitivity based on dry (P05) and wet (P95) references during the rice growing season for 5 classes in the Heilongjiang, Hunan, and Guangxi provinces.Due to the paddy rice flooded stage, the VH sensitivity is higher than that of other types.However, it still difficult to find a specific threshold at the province scale.

Figure 16 .
Figure 16.Violin plot of changes in VH sensitivity based on dry (P05) and wet (P95) references during the rice growing season for 5 classes in the Heilongjiang, Hunan, and Guangxi provinces.Due to the paddy rice flooded stage, the VH sensitivity is higher than that of other types.However, it still difficult to find a specific threshold at the province scale.

Table 2 .
The number of scenes from the Sentinel-1 C-band, Sentinel-2 MSI and Landsat 8 OLI that are used in this study in Heilongjiang, Hunan and Guangxi, respectively.

Table 2 .
The number of scenes from the Sentinel-1 C-band, Sentinel-2 MSI and Landsat 8 OLI that are used in this study in Heilongjiang, Hunan and Guangxi, respectively.

Table 3 .
The independent accuracy assessment of the rice extent product in Heilongjiang, Hunan, and Guangxi.

Table 3 .
The independent accuracy assessment of the rice extent product in Heilongjiang, Hunan, and Guangxi.