Complementarity of Two Rice Mapping Approaches : Characterizing Strata Mapped by Hypertemporal MODIS and Rice Paddy Identification Using Multitemporal SAR

Different rice crop information can be derived from different remote sensing sources to provide information for decision making and policies related to agricultural production and food security. The objective of this study is to generate complementary and comprehensive rice crop information from hypertemporal optical and multitemporal high-resolution SAR imagery. We demonstrate the use of MODIS data for rice-based system characterization and X-band SAR data from TerraSAR-X and CosmoSkyMed for the identification and detailed mapping of rice areas and flooding/transplanting dates. MODIS was classified using ISODATA to generate cropping calendar, cropping intensity, cropping pattern and rice ecosystem information. Season and location specific thresholds from field observations were used to generate detailed maps of rice areas and flooding/transplanting dates from the SAR data. Error matrices were used for the accuracy assessment of the MODIS-derived rice characteristics map and the SAR-derived detailed rice area map, while Root Mean Square Error (RMSE) and linear correlation were used to assess the TSX-derived flooding/transplanting dates. Results showed that multitemporal high spatial resolution SAR OPEN ACCESS Remote Sens. 2014, 6 12790 data is effective for mapping rice areas and flooding/transplanting dates with an overall accuracy of 90% and a kappa of 0.72 and that hypertemporal moderate-resolution optical imagery is effective for the basic characterization of rice areas with an overall accuracy that ranged from 62% to 87% and a kappa of 0.52 to 0.72. This study has also provided the first assessment of the temporal variation in the backscatter of rice from CSK and TSX using large incidence angles covering all rice crop stages from pre-season until harvest. This complementarity in optical and SAR data can be further exploited in the near future with the increased availability of space-borne optical and SAR sensors. This new information can help improve the identification of rice areas.


Introduction
Rice is a staple food for more than half of the world's population and the most important crop in low-income and lower-middle-income countries [1,2].Its contribution as a staple food and source of income and employment is crucial to food security [3].In 2012, rice was cultivated on about 163 million hectares of land that produced approximately 719 million tons of rice with 90% of it originating from Asia [4].The timely provision of accurate spatial information on the rice crop supports policies and decisions on agricultural production and food security in the region [5].
Remote sensing can provide spatial and temporal information to characterize rice agro-ecological attributes [6], assess rice growth [7] and productivity [8].It is also used to differentiate cropping patterns [9], detect cropping calendar [6] and the dynamic phenological stages of rice [10,11].
Passive and active sensors have different spatial, temporal and spectral characteristics that acquire different types of information at various levels of detail on the ground.Medium-resolution optical images such as SPOT [12,13], Landsat TM [14][15][16][17][18][19] and ETM [20] have been used successfully for paddy field delineation.Flood damage assessment in rice areas and detection of changes in rice area extent, composition and field conditions due to crop rotation, natural vegetation transformation and natural disasters (i.e., floods or storm) are other applications of medium-resolution satellite optical images [21].However, some studies based on SPOT-XS [13], Landsat TM and MOS-1 [22] reported difficulties in mapping the extent of rice area because of small field sizes and land cover heterogeneity in the area.Persistent and heavy cloud cover conditions also pose a problem when using relatively low-temporal-resolution images such as SPOT or Landsat TM [23], especially during the rainy season for vegetation mapping [24].
Optical images with low spatial resolution such as MODIS, AVHRR and SPOT VEGETATION provide wide area coverage with a high temporal resolution.These satellite images are suitable for capturing the dynamic phenology of rice [11,25] and mapping and monitoring changes in cropping systems [26], cropping intensity [27], cropping patterns [28] and cropping calendars [6,29] by season using vegetation indices such as NDVI data.Several studies on the use of MODIS conducted in China [30,31] and South and Southeast Asia [32] showed that MODIS is effective for mapping, monitoring and characterizing rice-based cropping systems for vast areas.
Active sensors such as Synthetic Aperture Radar (SAR) use longer wavelengths (0.01 to 1 m) to acquire ground information.SAR penetrates clouds and obtains information regardless of weather and light conditions [8].However, unusual weather such as intense localized events (rainfall or dense cloud cover) has been shown to affect high-frequency bands such as the X-band (2.4 to 3.8 cm) [33].Despite this, SAR remains a reliable sensor for providing stable data quality and images at a chosen time [34].
Studies have shown the successful use of SAR data for rice mapping, monitoring and yield estimation [23,[35][36][37].SAR's sensitivity to rice texture variations and flooded conditions, and its capability to capture the unique temporal backscatter signature of rice when observed throughout the crop season, makes SAR suitable for mapping cultivated rice area and monitoring growth [38].Capturing the flooding stage in SAR images is essential for the correct identification of rice fields [39].SAR images acquired at the right time to capture flooding [39] and early and late growth stages [40] improve rice crop identification.The availability of high-spatial-resolution SAR images with multitemporal resolution of 11 to 16 days such as Terra SAR-X and Cosmo Skymed (CSK) can be used to capture small patches of rice paddies.
Research on how the complementarity of both optical and SAR images can be used for comprehensive rice mapping to provide complete, timely and accurate spatial and temporal information on the rice crop is important in the agricultural sector.The use of low-resolution images such as MODIS can differentiate major land cover types [21] and provide an overview of the general situation on the ground.At a spatial resolution of 250 m, MODIS can be used for the delineation of land cover complexes with little information on the actual extent of individual land cover types.On the other hand, the use of a sensor with spatial resolution that is equal to or smaller than the size of the target fields can generate maps of currently irrigated croplands with high accuracy [41].Rice maps derived from high-spatial-resolution SAR images, however, require multidate imagery to be acquired at the proper time for high accuracy.
The objective of this study is to generate complementary and comprehensive rice crop information from hypertemporal moderate-resolution optical imagery and multitemporal high-resolution SAR based on remote sensing data and field data in the Philippines.The study relies on MODIS data to stratify and characterize rice areas at a small scale and X-band SAR data from TSX and CSK to identify cropped rice fields at a large scale and to extract details on the actual dates of flooding/transplanting during two consecutive seasons.

Study Area and Data
The study area covers the provinces of Nueva Ecija and Pangasinan located in Central Luzon, Philippines.These provinces are two of the main rice-growing areas in the Philippines that cover 25,878 km 2 (Figure 1).Nueva Ecija is the highest rice-producing province, while Pangasinan ranks third, with the highest rainfed production in Luzon.
Nueva Ecija is characterized by a terrain that comprises mostly low alluvial plains in the west and in the southwest, and rolling uplands in the northeast.The topography of Pangasinan, on the other hand, is characterized by broad alluvial plains in the central part and mountainous regions in the northeast and in the west [42].
The study area has a distinct wet season (WS) and dry season (DS), with maximum rainfall expected to occur from July to August.In general, the first or WS crop is sown in June, planted in late June or in July and harvested in late September or in October.The second or DS crop is sown or planted between late December and early January and harvested in April.Most rice farmers practice direct-seeding during the DS and transplanting in the WS.
Land preparation (plowing, harrowing, leveling) in Nueva Ecija and Pangasinan generally takes 3-4 weeks before planting.The majority of rice fields are prepared under flooded conditions in which fields have standing water of about 2-3 cm for about 3-7 days until the soil is soft enough for tillage [43].This is followed by leveling to remove high or low spots and ensure even water distribution and better weed control.Rice fields are usually planted 1 or 2 days after leveling.
In terms of rice production, irrigated rice is cultivated in areas where a national or communal irrigation system (NIS/CIS) is in place.These are large to medium-sized schemes supported by the National Irrigation Administration (NIA).Some small privately owned schemes also exist in the area.In general, two rice crops are cultivated annually in irrigated areas.Only one rice crop is cultivated during the WS in the remaining rainfed lowland areas, followed by fallow or a non-rice crop during the DS.There are very few upland rice areas characterized by unbounded fields with no water management.
Nueva Ecija is dominated by irrigated rice areas, which comprise 78% (138,157 ha) of the total physical area; the remaining 22% are rainfed (38,387 ha).In Pangasinan, 55% (98,313 ha) of the rice areas are irrigated, while 45% (82,001 ha) are rainfed [44], where the water supply is from rain and through groundwater pumping.Hence, the dry season rice area in Pangasinan is more fragmented than the more contiguous rice areas under surface water irrigation in Nueva Ecija.

MODIS Surface Reflectance 8-Day Composite
MODIS is a sensor onboard the Terra satellite launched by NASA into space in 1999 [45].The MODIS collection 5 product MOD09Q1 (Table 1) was used as the hypertemporal optical remote sensing data source for stratification and rice characterization, which is available at 8-day composite with 250-m spatial resolution.It is one of the MOD09 products that estimate ground surface reflectance corrected for the effects of atmospheric gases and aerosols.Each pixel contains the best observation values within the 8-day period with low view angle and the minimum amounts of clouds, haze or cloud shadow and aerosol loading [46].Each image contains red (band 1: 620 to 670 nm) and NIR (band 2: 841 to 876 nm), which were used to calculate NDVI (Equation ( 1)) [47]: CSK is a constellation of four X-band SAR satellites that became fully operational in 2010 under the mandate of the Italian Space Agency (ASI).The CSK Huge Region mode provides standard products with a resolution of 100 meters [48].However, an experimental product with 30-m spatial resolution has been processed (Table 1) and provided for the purpose of this study.The experimental data were acquired for DS 2012 to 2013.The revisit period using the same satellite with the same observation angle is 16 days.
TerraSAR-X (TSX) is a German X-band satellite operated by the German Aerospace Center (DLR) that was launched in June 2007 and became fully operational in 2008.TSX ScanSAR images were acquired in WS 2013.The revisit period is 11 days and the spatial resolution is 18.5 m [49] (Table 1).

Stratification
Stratified random sampling was used for the selection of rice monitoring sites.The classes derived from a MODIS unsupervised classification were used as the strata for the random selection of monitoring sites (c.f.Section 2.2).A pre-test of the survey questionnaire (Figure S1) and the associated fieldwork activities (including GPS coordinate collection and site observations) and the travel time to go from one monitoring site to another were used as the basis to determine the number of samples that could be accomplished in one day.Budget, person-days, labor constraints and nearness to the road were additional considerations.Based on these, a sample size of 11 per class was used which was derived using the sample size equation (Equation ( 2)):  2).GPS coordinates were collected for each site.Interviews with farmers (Figure S1) were carried out to collect data on cropping systems, crop calendar (a sequential summary of periods of operations such as land preparation, planting and harvesting) [50] and rice ecosystem (irrigated or rainfed) [51].
Cropping system included sequential cropping information such as cropping intensity (the number of crops planted on the same field within a cropping year: single, double or triple) [6] and cropping pattern (temporal sequence of crops planted in an area for the WS and DS within a cropping year [52]: rice-rice, rice-other, other-rice, rice-fallow, fallow-rice).Field monitoring was carried out to determine the actual condition of rice fields at the time of every SAR acquisition.Separate data sets during the DS (33 fields) and WS (36 fields) were monitored (Table 3).The selected rice fields were observed and the field conditions such as land preparation status (plowed, harrowed, leveled, flooded, saturated), crop establishment, growth stages and weather condition at the time of the field visit were recorded.Actual dates of flooding and transplanting were acquired through interviews with farmers for the WS 2013.

MODIS-Based Rice System Stratification and Characterization
Figure 2 shows the integrated flow chart for rice stratification and characterization using MODIS data and rice crop identification using SAR data.Hypertemporal MODIS NDVI images were processed to generate strata that reflect differences in practiced cropping systems.
MODIS 8-day composite images from 2006 to 2011 were used for deriving 276 NDVI images.Linear stretching was applied to generate the digital NDVI numbers (DN), where −1, the minimum NDVI value, was assigned to 0 and the maximum NDVI value of 1 to 255.The upper envelope of the NDVI time series was derived using the adaptive Savitzky-Golay filter in TIMESAT [53][54][55] (Figure S2) and the resulting time series was classified using the ISODATA [56] clustering algorithm through unsupervised classification [9,[57][58][59].
Each ISODATA run was evaluated using the divergence distance measure [60] that assesses the clustering signature separability.The highest positive deviation from the trend in average divergence indicates a solution in number of classes that can be extracted from the time series [9,59].The extraction of 115 classes provided an optimal stratification for the NDVI time series.As described, the temporal signatures of the 115 NDVI classes (Figure S3) were evaluated for merging the likely non-rice class into one.
Differentiation of the rice and non-rice strata based on their temporal NDVI profiles was carried out using the existing 1980s rice area map, visual interpretation of Google Earth imagery, experts' knowledge and opinion exchanges among authors.The field data were randomly divided into two parts.One part (144 points) was used to define the classes, that is, to assign the corresponding rice-based cropping system characteristics (collected on the ground and from farmers' interviews) to each of the eight NDVI classes.(1) Strip mosaicking and multilooking-single frames in slant range geometry with the same orbit were mosaicked along their azimuth to facilitate data processing and handling.Multilooking was carried out to improve SAR image quality through a reduction in speckle and to obtain approximately square pixels [62].Averaging of range and resolution cells generates the multilooked images [63].

SAR-Based Rice Crop Identification and Planting Period Mapping
(2) DEM-based orbital correction-SRTM 90 m Digital Elevation Model (DEM) tiles were used for the orbital correction.Errors of the azimuth start time and/or slant range distance were corrected on the basis of a reference DEM.
(3) Co-registration-strip mosaics covering the same area with the same geometry and mode were co-registered by gross shift estimation based on orbital data.A set of sub-windows based on reference image and images to be used for co-registration was then automatically identified.The shift between pixels, including elevation, was calculated through cross-correlation.Polynomial function was used to calculate the shifts to be applied in azimuth and range direction [61].
(4) De Grandi time-series filtering-a balance of differences in reflectivity between images at different times was achieved with the use of an optimum weighting filter [64].Multitemporal filtering is based on the assumption that SAR geometry is the same for all acquisitions.The reflectivity can change because of dielectric and geometrical properties but not because of a different position of the resolution element with respect to the radar [61].
(5) Terrain geocoding, radiometric calibration and normalization-conversion of backscatter elements into slant range image coordinates was carried out using a DEM as a backward solution.Range-Doppler equations [65] were used in the transformation of two-dimensional coordinates of the slant range image to three-dimensional object coordinates in a cartographic reference system.Geometric and radiometric calibration of the backscatter values is necessary for inter-comparison of radar images acquired at different times with different sensors and/or different viewing geometries [63].Radiometric calibration was performed using a radar equation that takes into consideration the scattering area, antenna gain patterns and range spread loss.The backscatter coefficient (σ°) was normalized to compensate for the range dependency using the cosine law of the incidence angle [61].
(6) Anisotropic Non-Linear Diffusion (ANLD) filtering-the ANLD filter performs strong smoothing in homogeneous areas while preserving signal variations coming from neighboring areas [61] and linear structures (e.g., roads, rivers and field edges).A diffusion equation was used wherein the diffusion coefficient is a function of image positions and assumes a tensor value [66].
(7) Removal of cloud-related effects from localized intense weather events-anomalous peaks or troughs caused by localized intense events were identified through the analysis of temporal σ°, which was corrected using an interpolator.A priori information on the cropping calendar and weather conditions at the time of image acquisition is necessary for correct interpretation of these events [33,61].
These processing steps were carried out in an almost fully automatic way to allow the analysis and interpretation of the multitemporal backscatter signature.

Rice Detection Algorithm and Threshold Selection
After the SAR data processing, specific thresholds (Table 4) were set in MAPscape-RICE to identify and map rice as well as the flooding/transplanting dates for each season.Flooded soils, which are a unique characteristic of paddy rice fields [7,30], are characterized by low SAR backscatter that makes them easily separable from other areas with no surface water and certain ground roughness.Stable water has a low average backscatter (i.e., a low mean value measured across the whole multitemporal series) compared with rice.Built-up areas, on the other hand, show higher average backscatter than rice [34].In general, rice can be discriminated from other crops that are not inundated during the planting and vegetative stage if the multitemporal acquisition series is properly synchronized with the cropping calendar.The minimum rice cycle duration is 90-100 days [67] increasing to 140-150 in temperate latitudes [68].The majority of rice farmers plant a 110-120-day rice variety in the Philippines.Figure 3 shows examples of temporal backscatter signatures for different land covers.The algorithm used to identify the rice crop using multitemporal images is very similar to the algorithm that is fully documented in [61] and the parameters (P) we used for rice detection are shown in Table 4.We first excluded non-rice areas such as permanent water bodies (P5) where there is a stable low backscatter and built-up areas (P6) where there is a stable high backscatter (Figure 3).Areas that exhibit non-permanent, longer duration flooding (P7) such as fishponds, irrigation tanks and seasonal wetlands and other seasonal field crops that show larger variation (P4) in backscatter from the growth in biomass were also masked out.The remaining unmasked pixels undergo a stepwise process to look for agronomic flooding where a drop, that is, the lowest value in the temporal signature of rice is, observed (P1).This marks the start of the rice growing cycle (SoS) and is also identified as the flooding/transplanting date.Once agronomic flooding is detected, we look for a rapid increase in biomass that confirms that rice is grown.The parameters used to detect rice growth are PoS (P2), span of SoS to PoS (P3), minimum growth (P4) and minimum rice cycle duration (P8).If the flooding period was not captured in any of the image acquisitions, the same parameters in detecting rice growth are used [61].The non-rice area exclusion and rice growth detection through a rapid increase in biomass are also similar to the criteria used for rice paddy mapping using MODIS time-series data [32].
Selection of thresholds was guided by the backscatter coefficient extracted from the monitored sites.For every image acquisition, field visits were carried out at the selected sites to link the backscatter values with the corresponding land preparation condition and growth stages.From this, simple statistics such as the minimum, maximum, mean and range were calculated to set the ceiling for the threshold selection.Table 4 describes the threshold parameters used for the identification of rice areas and flooding/transplanting dates.
The rice identification using thresholds generated rice maps with 15-m pixel size for CSK and 10 m for TSX (these values correspond to the pixel sampling of the original SAR products).The pixel size follows the Nyquist sampling theorem, wherein a sampling rate higher than the real pixel resolution (i.e., half of the resolution cell) [69] is required for "perfect" interpolation fidelity.

Accuracy Assessment
Accuracy assessment of the MODIS-derived rice map was based on the kappa coefficient and confusion matrix (Tables S1-S8) showing the correctly classified pixels with reference to ground-truth points [70].The rice map generated from MODIS was evaluated using 109 points (the rest of the points were used for defining the classes derived from ISODATA clustering) with ground information indicating whether an area is rice or non-rice.Assessment of rice cropping intensity required ground information on the number of times an area is planted with rice in a year for the WS and DS (single or double rice), while cropping intensity is compared to information on the number of crops planted in a year regardless of whether it is rice or another crop.Cropping pattern accuracy assessment used the sequential pattern indicating the type of crops planted in a year (rice-rice, rice-fallow, rice-other, fallow-rice).The crop calendar for WS and DS was assessed separately using ground information on the general planting dates.
The spatially detailed rice area map generated from SAR was evaluated using 240 points of rice/non-rice ground information that covers the TSX footprint (Table 3).Root Mean Square Error (RMSE) was used to determine the difference between the actual flooding/transplanting dates (from farmers' interviews and ground observations) and the dates derived from TSX.An RMSE of less than 11 days indicates good detection of planting dates using TSX ScanSAR multitemporal images.The relationship between the actual planting dates and the TSX-derived dates was assessed using the coefficient of determination (R 2 ).DS flooding/transplanting dates derived from CSK ScanSAR were not assessed for their accuracy due to the lack of ground information on the actual planting and sowing dates.

MODIS-Based Rice Stratification and Characterization
A rice characteristics map was generated using ISODATA unsupervised classification of MODIS-derived NDVI images.From the preliminary result of 115 classes, all non-rice classes were merged into one, and the remaining possible rice classes with similar NDVI temporal profiles (Figure S3) were merged based on the visual assessment of signature similarity.This resulted in a map with 13 rice classes and one non-rice class.Field data on the characteristics of different rice cropping systems were used to describe and assign cropping systems to the MODIS-derived NDVI strata or classes.The 13 rice classes were merged into nine classes (Figure 4) based on similarity in cropping calendars.The final stratified rice map contains the following information on rice cropping systems: cropping calendar, cropping intensity, cropping pattern and rice ecosystem (irrigated or rainfed).
Accuracy assessment of the rice-based cropping system characteristics map using the ground data showed an overall accuracy of 87.2% for rice or non-rice and 82.7% for cropping intensity (Table 5).The confusion matrix is in the supplementary information Tables S1-S8.In the case of rice ecosystem discrimination (74.1%), the result shows that, inside the boundary of NIA irrigation divisions, some rice fields are not receiving irrigation water and they are thus identified by farmers during the interview as rainfed.As we mapped more details such as rice cropping intensity (77.8%), cropping pattern (77.8%) and cropping calendar (WS: 75.3% and DS: 61.7%), the overall accuracy became lower.In summary, the results showed an overall accuracy of 61.7% to 87.2% with a kappa coefficient of 0.52 to 0.72 (Table 4).Higher accuracies were obtained for the most basic information layers, such as rice/non-rice or cropping intensity, and lower accuracies were obtained for the more thematically detailed information layers such as the cropping calendar.

Temporal Rice Backscatter Signature from TSX and CSK ScanSAR and Threshold Selection
Monitoring of the rice areas during the DS 2012-2013 and WS 2013 for every satellite image acquisition facilitated establishment of the rice temporal backscatter signature from pre-season to harvesting for CSK and TSX ScanSAR, respectively.Figure 5 shows the range of backscatter values extracted from CSK and TSX data for the monitoring sites, grouped by the observed land preparation and rice crop growth stage.
Tables 6 and 7 contain the thresholds used for rice identification and flooding/transplanting dates.The rice temporal backscatter signature was used to guide the threshold selection.Simple statistics such as the mean, minimum, maximum and range were used in setting the ceiling for the threshold selection.A description of the parameters is given in Table 4. for each observed land preparation and growth stage at the monitored sites.The black dots represent the outliers, the thick horizontal black line in the middle of the box is the median, the upper half of the box is the 25th percentile and the lower half is the 75th percentile, and the extent of dashed lines represent the minimum and maximum.

SAR-Based Rice Cropped Area and Flooding/Transplanting Dates Map (TSX and CSK ScanSAR)
The temporal backscatter signatures observed in Figure 4 were used to select specific thresholds that were applied to the multitemporal SAR data to generate rice area maps and maps of planting dates for the DS 2012-2013 and WS 2013.The combination of these consecutive seasonal area maps resulted in a physical rice area map (Figure 6) containing information on areas planted with rice across seasons.The overall accuracy was 90.4% with a kappa coefficient of 0.72 (Table 5).
The spatially detailed flooding/transplanting dates for the DS and WS generated from SAR are shown in Figures 7 and 8     Accuracy of the flooding/transplanting dates was assessed by linear correlation.The result showed strong correlation (R 2 = 0.87) between the TSX-derived dates and the ground truth dates (Figure 9).Root Mean Square error (RMSE) was also calculated to determine the difference between the TSX-derived and the actual flooding/planting dates.The RMSE is 9 days (Table 8), which is less than the 11-day revisit period of TSX ScanSAR.

Comparison of MODIS and SAR Planting Dates
Figure 10 shows the comparison of the usual/general planting period acquired from MODIS and the planting dates from ScanSAR for DS 2012 and WS 2013.The DS is characterized by a much shorter time period, reflecting the narrow water release dates in the irrigation systems.The WS shows a much larger range of planting dates, albeit with a defined peak planting window.A shift in the planting calendar can be observed where the planting period from CSK and TSX exhibited a late (samples A and F) and early (samples D and E) start of planting with reference to the general planting period from MODIS.Sample D also shows that planting dates in August from TSX were not captured as there was no image acquisition during that month.Sample B, on the other hand, showed no difference for MODIS and CSK.Sample C showed a different result.In MODIS, this rice cropping system does not plant rice in the WS; however, TSX shows that some pixels were detected as rice in WS 2013 but with no obvious peak planting period.This could be attributed to the higher spatial resolution of SAR that was able to pick up rice in these relatively small areas.

MODIS-Based Rice Stratification and Characterization
The MODIS-derived rice characteristics map contains information on the distribution of the predominant cropping systems being followed in the study area.MODIS was found to be effective in capturing the large homogeneous rice areas (B, C and E in Figure 3).However, rice areas with small field sizes and variability in the cropping systems practiced such as the cropping calendar and crops planted mostly located in the rainfed areas tend to be generalized within the strata when using MODIS.The identification of cropping systems in the smaller and dispersed rice classes (i.e., D, G and H from Figure 3) in the rainfed areas was found to be challenging.Overall, we note that the ability of hypertemporal MODIS to map the predominant rice characteristics is related to the spatial homogeneity and spatial fragmentation of the observed area.

Nov-Dec
Dec-Jan

Fallow during WS from MODIS Aug
May-Jun Jun-Jul

Rice SAR Temporal Backscatter Signature and SAR-Based Rice Cropped Area and Flooding/Transplanting Dates Map
SAR sensors (TSX and CSK) were found to be effective in mapping the rice areas with an overall accuracy of 90.4%.Newly planted or transplanted rice showed an increase in backscatter coefficient from both CSK and TSX (Figure 5) at an incidence angle of 40° and 45°, respectively.This result confirms the findings of [7] that high-frequency bands such as Ka, Ku and X bands can detect thin rice seedlings after transplanting at large incidence angles [7].During the initial part of the growing season, rice shows an increase in backscatter because of a progressively higher surface roughness (i.e., from transplanting until tillering to stem elongation stage).As the rice canopy closes to cover the stems and spaces between rows and the surface becomes more homogeneous, the backscatter decreases because of a progressively lower surface roughness.High-frequency bands are sensitive to early rice growth stages (with LAI < 1.5) within one month after transplanting and lose sensitivity after that [7].Backscatter of rice from X-band saturates at an earlier growth stage [34,71].This research has identified the peak to occur within the tillering to stem elongation stage.This is then followed by a decrease in backscatter coefficient that starts from the panicle initiation stage.An increase in backscatter was observed to occur again during the maturity stage, which is also in agreement with the findings of [7] that X-band is sensitive to ripening grains in the larger incidence angle.Harvested areas during the WS show higher backscatter than in the DS (Figure 5).This can be attributed to the moist condition related to the dielectric properties of stubbles during the rainy season (WS), whereas, during the DS, very dry stubbles and soils exist in the harvested areas.
Multitemporal X-band SAR data clearly identifies newly transplanted rice seedlings and hence has the potential to detect transplanting dates [34].The analysis of the relationship between TSX-derived dates and the actual flooding/transplanting dates showed a very strong correlation with an R 2 of 0.87 (Table 7).CSK-derived flooding/transplanting dates were not assessed in the same way because of a lack of information on the actual flooding/transplanting dates.
Factors affecting the accuracy of the SAR-generated rice maps include weather conditions and the timing of image acquisition.Several tropical storms passed over the study area during the WS acquisition period for TSX ScanSAR, which affected the images.Examples of the effect of these storms are shown in the Supplementary Information section (Figure S4).We assume that thick cloud cover in the image appears as black patches that decreased the intensity of backscatter from the ground that was not able to penetrate the clouds to reach the SAR sensor.The intense localized rainfall, on the other hand, is assumed to be the bright areas in the image that were probably produced because of the increase in water in the air that acted as multiple scatterers, also increasing the dielectric constant, resulting in an increase in backscatter.Further ground data are needed to validate these observations, which were not available in this study since the intense rainfall precluded any field activity at the time and location of these events.Longer microwave bands (e.g., C-band) for rice mapping, especially during the wet season, should be considered when mapping rice areas in a tropical country such as the Philippines, which experiences around 20 severe tropical storms every year [72].Proper planning and timing of the image acquisitions are also essential for accurate rice area mapping.In our methodology, the rice area will be strongly underestimated if the SAR time series does not adequately capture the flooding/transplanting dates.
In our study, one image was not acquired in August, which is toward the end of the transplanting period, and it is likely that we underestimated the area of the rice crop that was established around that time.

Comparison of General Planting Dates from MODIS with SAR-Derived Flooding/Transplanting Dates
Our comparison of the general planting period from MODIS with the SAR-derived flooding/planting dates showed consistency in the periods of planting for sample B and a slight early shift of planting for samples A, E and F. The use of MODIS time series (2006-2011) for capturing the general planting period was shown to be a valid approach since the majority of the rice areas detected by SAR in 2012-2013 were planted in time periods that were consistent with the general planting dates derived from MODIS (Figure 10).

Conclusions
This study showed how hypertemporal MODIS and multitemporal, high-spatial-resolution SAR images can be used to derive complementary information from rice areas.We used hypertemporal MODIS data to generate a rice characteristics map of the dominant rice-based cropping systems containing information on the general cropping calendar, cropping pattern, cropping intensity and rice ecosystem, while X-band SAR was used for the spatially detailed and accurate mapping of rice areas and flooding/transplanting dates.
The study also provided the first assessment of the temporal variation in backscatter in X-band SAR from CSK and TSX using large incidence angles (40° and 45°, respectively) for rice covering stages from pre-season and land preparation until harvest.CSK and TSX with large incidence angles can detect newly transplanted seedlings.An increase in backscatter was observed for both CSK and TSX from flooding to transplanting.The peak in backscatter coefficient was observed to occur from tillering to stem elongation, which is followed by a decrease in backscatter due to lower surface roughness as the canopy closes and covers the stem and spaces between rows.An increase in backscatter was observed again at the maturity stage, with the ripening grains as the major scatterers with large incidence angles.In comparison to C-band SAR such as ERS-1, flooded rice fields have increasing backscatter until they reach the reproductive phase, which can be attributed to the wave-vegetation-water interaction [68].C-band SAR backscatter reaches the maximum at early ripening phase, which is attributed to the vertical plant structures and water surfaces as scatterers, while a slight decrease in backscatter occurs during the late ripening phase due to a reduction in the capability to penetrate the soil or water [73].
High accuracy (90.4%) was obtained from the spatially detailed rice area map derived from SAR.The flooding/transplanting dates from TSX ScanSAR had a very strong correlation with the actual planting dates, showing its effectiveness in detecting the flooding/transplanting period.Accuracy of the information that describes the characteristics of the cropping systems within the rice areas obtained from MODIS ranged from 61.7% to 87.2%.The most basic information layers (e.g., rice/non-rice) had the highest overall accuracy, while the more thematically detailed information layers such as the cropping calendars had the lowest overall accuracy.MODIS provides general information that leads to a better understanding of rice areas, which gives insight into general farmers' management practices within each rice cropping system.
The results suggest that the planting period derived from MODIS can be used as a guide in planning SAR image acquisition [41] even over large extents in addition to secondary information on irrigation water release dates.MODIS can provide a general planting period while SAR can provide spatially detailed planting dates.Though X-band SAR is effective in mapping the temporal variability in planting dates and in detailed mapping of rice areas, the key to an accurate mapping is to acquire SAR images at the right time and with a sufficiently high repetition frequency, especially during the critical phases of the cropping season.This includes the flooding period, which is needed for the correct identification of rice areas [39], and the tillering to stem elongation stage, when the maximum backscatter is achieved, which confirms that rice is planted and grown.Missing an image during these critical periods could lead to unmapped rice area as experienced in this study.
This study emphasizes the complementarity of optical and SAR sensors to acquire comprehensive rice information.Together, the use of optical and SAR sensors generated accurate, comprehensive and necessary information to help aid decision making in the identification of rice areas for intensification, and areas for the development of irrigation as one of the necessary steps in dealing with food security issues.Despite some limitation in the SAR acquisitions, this study suggests that current and future optical and SAR systems such as Sentinel 1, 2, 3, Proba-V, RISAT-1, Landsat 8, MODIS, CSK and TerraSAR-X can be used complementarily to provide valuable spatial, thematic and temporal information on the rice crop in Asia.

Figure 1 .
Figure 1.Pangasinan and Nueva Ecija showing rice area in the 1980s (green).The brown and blue boxes delineate coverage of CSK and TSX SAR products, respectively.MODIS tile h29v07 completely covers the area.The yellow line shows the boundary of irrigated areas as reported by NIA (2011) and the sites for ground data collection are in black dots.

(8 20 4 .
) /14 number of points in one day number of working days sample size number of classes from MODIS ISODATAunsupervised classification Field Data Collection During the 2012-2013 DS, 253 sites were visited, that is, 194 sites in rice areas and 59 in non-rice areas (Table

Figure 2 .
Figure 2. Flow chart of data processing method.The asterisk (*) indicates the use of MODIS in determining the dates for acquisition of ScanSAR images.

2. 3 . 1 .
Basic Processing of Multitemporal SAR Images SAR images were acquired over the study area for the DS (2012 to 2013) and WS (2013) to identify rice areas and determine flooding/transplanting dates based on time-series analysis of the backscatter coefficient (dB).The set of CSK and TSX ScanSAR images in Single Look Complex (SLC) format was processed separately using MAPscape-RICE, a semi-automated software developed by sarmap.The SLC data were transformed into a terrain geocoded backscatter coefficient by following these processing steps[33,61]:

Figure 3 .
Figure 3. Example of a temporal backscatter signature from different land covers derived from CSK (left) and TSX (right).

Figure 4 .
Figure 4. MODIS-derived rice cropping systems and characteristics map.

Figure 5 .
Figure 5. Boxplots of the backscatter coefficient derived from CSK (left) and TSX (right)for each observed land preparation and growth stage at the monitored sites.The black dots represent the outliers, the thick horizontal black line in the middle of the box is the median, the upper half of the box is the 25th percentile and the lower half is the 75th percentile, and the extent of dashed lines represent the minimum and maximum.
. The flooding/transplanting dates that were detected during satellite acquisition for the DS were from 6 December 2012 to 8 February 2013.For the WS 2013, the dates were from 25 May to 30 July 2013.

Figure 8 .
Figure 8. Flooding/transplanting dates within the rice areas during WS 2013 generated from TSX ScanSAR multitemporal images.

Figure 9 .
Figure 9. Correlation of TSX-derived planting dates with the ground data.

Figure 10 .
Figure 10.Examples of planting period from different rice cropping systems acquired from MODIS (blue box) and SAR (gray columns) for the DS and WS.

Table 1 .
Characteristics of the series of satellite images acquired for rice mapping.

Table 2 .
Number of rice fields and non-rice areas visited for data collection.

Table 3 .
Number of rice fields used for the accuracy assessment of rice-based cropping system characteristics and rice area delineation.

Table 4 .
Threshold parameters for rice identification and detection of flooding/transplanting dates.

Table 5 .
Accuracy assessment of the cropping system information derived from MODIS for the whole Nueva Ecija and Pangasinan provinces.

Table 6 .
Thresholds selected for rice identification and flooding/transplanting dates using TSX multitemporal images.

Table 7 .
Thresholds selected for rice identification and flooding/transplanting dates using CSK multitemporal images.

Table 8 .
Correlation of TSX-derived flooding/planting dates with the ground truth flooding/planting dates.