EXPLOITING SENTINEL-1 SAR TIME SERIES TO DETECT GRASSLANDS IN THE NORTHERN BRAZILIAN AMAZON

Recent advances in cloud-computing technologies and remote sensing data availability foster the development of studies based on the analysis of optical and SAR imagery time series. In this paper, we assess the potential of Sentinel-1 imagery time series for grassland detection in the northern Brazilian Amazon. We used the Google Earth Engine cloud-computing platform as an alternative to obtain and analyse Sentinel-1 imagery, acquired from 2017 to 2018 over the region of Mojuí dos Campos/PA, Brazil. We extracted several temporal metrics from the imagery time series and used the Random Forest algorithm to perform the classification. In addition, we analysed the time series considering different channels, including the VV and VH polarizations, both separately and in combination, and the CR, RGI and NL indices. We could efficiently discriminate areas of grasslands from forest and agricultural crops using either VH time features or features extracted from the combination of both VV and VH polarizations. The classification map that resulted from the combination of VV and VH data presented the highest accuracy, with an overall accuracy of 95.33% and a 0.93 kappa index. Despite simple, the approach adopted in this paper showed potential to differ grasslands from areas of agriculture and forest in the northern Brazilian Amazon.


INTRODUCTION
The dynamics of deforestation observed in the Brazilian Amazon are complex and, for that reason, have motivated several studies that investigate land use and land cover (LULC) changes. For instance, there are studies that investigate factors that contribute to the Amazon deforestation (Barona et al., 2010), the LULC change dynamics of deforested areas (Carreiras et al., 2014) and, its environmental implications (Farias et al., 2018). Barona et al. (2010) indicates agricultural crops and grasslands as the main land cover types that occupy the deforested areas of the Brazilian Amazon. In fact, the authors suggest that the pasture expansion might contribute more to the deforestation of the Brazilian Amazon than agricultural activities, such as the soy industry. Some reports also support this statement, showing that grasslands occupy most of the Brazilian Amazon deforested areas (EMBRAPA and INPE, 2011;Veiga et al., 2004;Zero Deforestation Working Group -ZDWG., 2017).
The Brazilian Agricultural Research Corporation (EMBRAPA) and the Brazilian National Institute for Space Research (INPE) released an executive summary in 2011 showing that over sixty percent of the areas deforested until 2008, in the Brazilian Amazon, shifted to grasslands. In addition, ZDWG (2017) also presented, in 2017, that sixty-five percent of the Brazilian Amazon deforested areas are used for low-efficiency pastures. These aspects highlight the importance of grassland sustainability to agricultural activities in the Amazon. Grasslands occupy a significant portion of the agricultural area of the Brazilian Amazon, contributing to agricultural and livestock systems (Navegantes-Alves et al., 2012;Taravat et al., 2019). These aspects emphasize the importance of the detection and monitoring of grasslands. * Corresponding author In this context, remote sensing techniques are commonly used, due to their potential to capture data useful to describe LULC changes. However, LULC mapping might be challenging because of the limitations associated with the acquisition and processing of satellite imagery time series (TS), as the frequent cloud cover observed during the rainy season in tropical regions reduces the availability of cloud-free imagery from optical sensors. In this sense, the use of Synthetic Aperture Radar (SAR) technologies represents an alternative to the acquisition of satellite imagery over regions that present high cloud frequency. SAR sensors provide the opportunity to acquire consistent remote sensing data at day and night, without being as much affected by cloud and weather conditions as optical sensors.
In addition, the analysis of large volume data require high computing power and storage capability (Carrasco et al., 2019), which used to represent a common issue. However, the advent of cloud-computing infrastructures supports data integration across different data centres in scalable and flexible platforms with significant computing power. Carrasco et al. (2019) highlight that platforms such as Google Earth Engine (GEE) (Gorelick et al., 2017) foster the development of new methods for LULC mapping, once it allows applying complex image processing algorithms to large amounts of data. In the context of LULC change detection, the analysis of dense satellite imagery TS might contribute to the identification of environmental disturbances that could pass unnoticed if used short imagery TS. In tropical regions, it is possible to observe this effect due to the rapid development and recovery of the vegetation (McDowell et al., 2015;Shimizu et al., 2019).
The use of cloud-computing platforms has attracted great attention of the remote sensing community in the last years, particularly for the detection of LULC changes, based on imagery TS analysis (Carrasco et al., 2019;Lee et al., 2018;Mardani et al., 2019;Sidhu et al., 2018). In this scenario, the use of these platforms to access and process SAR TS has also provided support to studies related to the monitoring of tropical forests (Chen et al., 2018(Chen et al., , 2017Shimizu et al., 2019). These technologies show potential to contribute to several remote sensing researches that require accessing and processing SAR imagery TS.
In this study, we address the problem of grassland detection in the Brazilian Amazon, based on the analysis of SAR imagery TS. We used the GEE platform as an alternative to obtain and analyse Sentinel-1 imagery, acquired from 2017 to 2018 over the region of Mojuí dos Campos, Pará (PA), Brazil. The algorithm extracts several multitemporal features from the TS and uses the Random Forest (RF) classifier to identify the patterns that characterize grasslands. In the remote sensing-related literature, few studies focus on the use of SAR imagery TS to detect grasslands (Taravat et al., 2019;Whelen and Siqueira, 2018). This is one of the first studies directed to the detection of grasslands in the Brazilian Amazon based SAR TS analysis.

Study Site
The study area of the proposed method is mostly located at Mojuí dos Campos municipality, including part of the Belterra municipality, both at Pará, Brazil (Fig. 1). The study site comprises about 771 km², being mostly composed of agricultural crops, pasture areas and native forest. The Köppen climate classification (Köppen and Geiger, 1928) defines the study site as a region characterized by a tropical monsoon climate (class 'Am'). The annual average temperature in the region tends to vary from 27°C to 30°C, with a monthly precipitation of about 184mm during the dry season (April to October), and approximately 239mm during the wet season (November to March) (Alvares et al., 2013).

Sentinel-1 Data
The Sentinel missions, developed by the European Space Agency (ESA) under the Copernicus Programme, represent a satellite family composed of different sensors, including multi-spectral and radar imaging instruments. Each Sentinel mission is composed of a constellation of two satellites in order to achieve the ESA revisit and coverage requirements (ESA, 2017).
Sentinel-1 represents a polar-orbiting radar imaging mission that provides all-weather and day-and-night data, suitable for both land and ocean services (EOS, 2014).
The Sentinel-1 satellites, launched on April 2014 (Sentinel-1A) and April 2016 (Sentinel-1B), include sensors that operate at Cband wavelength, with different acquisition modes, including the Interferometric Wide Swath (IW). This acquisition mode is suitable for land cover applications with a 5m x 20m resolution in range and azimuth, respectively. In this mode, the sensor acquires data in dual-polarization, which can be HH+HV or VV+VH, with a wide swath of 250 km and a high revisit time (12 days, considering only one satellite) (ESA, 2013) The Sentinel-1 imagery TS analysed in this study is composed of 31 IW scenes, obtained with VV and VH polarizations, from September 22 nd , 2017, to September 17 th , 2018. In addition, we only considered S1A descending orbit data, in order to standardize the analysed time series. Some pre-processing steps were also applied by the GEE platform in order to derive the backscatter coefficient ( ) in each pixel and to provide direct access to orthorectified images. We can describe as the normalized measure for the intensity of each pixel. The pre-processing steps include border and thermal noise removal, radiometric calibration and terrain correction (GEE, 2018). In the terrain correction procedure, GEE converts data from ground range geometry to using the 30m-resolution SRTM Digital Elevation Model (DEM) (GEE, 2018). The fast and efficient access to pre-processed Sentinel-1 SAR imagery, provided by GEE, promotes the development of different remote sensing studies. Although the backscatter coefficient offers a convenient way to use and process SAR data, it does not convey phase information, which is lost during computation. Figure 2 presents an overview of the methodology. First, we applied a multitemporal filter developed by Quegan and Yu (2001) followed by a 7x7 Frost filter (Frost et al., 1982) to each image that composes the pre-processed imagery TS. The multitemporal speckle filter assumes that the images of the time series are sequentially uncorrelated (Maghsoudi et al., 2012). This filter, which results in speckle-reduced images, uses the input intensity data and an estimate of the local mean backscattering coefficient as entrance data. The multitemporal speckle filter follows the following formula:

Methodology
where , = radar intensity of the output image k at pixel , 〈 , 〉 = estimate of the local mean backscattering coefficient , = input intensity data M = number of images 〈 , 〉 = spatial average over the N pixels in a local window that surrounds the current , position in image The Frost filter allows spatially reducing multiplicative noise from images, while the multitemporal filter consider and minimizes the speckle influence on each pixel of imagery time series. The Frost filter is a convolutional filter that works with an n-by-n moving kernel that replaces the central pixel by the weighted sum of the neighbouring pixels in kernel (Wang et al., 2012). The Frost filter follows the following formula: After these pre-processing procedures, we computed three different radar indices, taking into account the information provided by both polarizations. These indices are: the VH/VV polarization ratio, defined here as cross ratio (CR); a normalization ratio (NL) (Lu et al., 2011); and, the Radar Gap Index (RGI) and are given by: At this point, we can separate the imagery dataset into five time series, accordingly to the attributes (VV, VH, CR, NL and RGI).
In this study, we analysed separately each one of these time series. In order to do so, we extracted several temporal metrics, including minimum, mean, maximum, standard deviation, sum, amplitude -defined as the difference between the maximum and minimum value -and coefficient of variation. In addition, we also performed the classification using metrics extracted from the combination of both VV and VH time series. The classification method used to classify these temporal metrics was the Random Forest (RF) classifier. This classification method generates a large number of decision trees, based on the use of different subsets of training samples. This algorithm, widely used in the remote sensing community over the last decades (Belgiu and Drăgu, 2016;Diniz, 2019;Sica et al., 2019), is reported to provide high accuracy classification results (Breiman, 2001). In this study, we set the classifier to grow 300 trees (ntree) and to use the square root of the number of the predictor variables for splitting at each tree node (mtry).
The input data required by the RF algorithm include not only the temporal metrics but also training samples. We used the MapBiomas LULC mappings (MapBiomas, 2018) related to the last three years (2016, 2017 and 2018) in order to obtain a sample dataset composed only of invariant samples. We adopted a stratified sampling strategy to define 500 random points that represent the pattern associated with the classes Agriculture, Forest and Grassland. For that reason, we considered as sampling areas only regions of agriculture, forest and pasture larger than 50 hectares that presented the same classification in these three LULC mappings. We split the sample collection into two datasets, considering 60% of the samples for training and 40% for testing the classification model.
In order to avoid noises in the classification results, the RF classifier was integrated with the majority filter operation, considering a 3x3 kernel. This contextual filter can efficiently reduce classification noises commonly observed in products obtained by pixel-based classification approaches. This method employs a moving window over the entire classification map. It reclassifies the central pixel of the moving window to the class attributed to the majority of the pixels contained in that position of the window. In cases where there is not a main class in the window, the algorithm maintains the class of the central pixel.
To validate the classification maps generated by the RF classifier using the temporal metrics, we used the testing dataset composed of invariant samples. We analysed statistics provided by the confusion matrix, including the overall accuracy (OA) and the kappa index. In addition, we compared the statistics obtained for the classification maps resulting from the use of the VV, VH, CR, RGI, NL and VV+VH time series. We used the Z-test as an alternative to evaluate whether there are significant differences between these products. Figure 3 illustrates the invariant samples of agriculture, forest and grassland (500 per class), considering a false-colour composite image that uses the metrics extracted from the VH time series. In this figure, we consider the coefficient of variation, minimum and amplitude images to compose the Red, Green and Blue (RGB) channels, respectively. Note that this false-colour composite emphasizes the differences between the analysed classes. In this case, green pixels put in evidence areas of forest, with dense vegetation cover, while pink and purple pixels indicate areas of agriculture and grasslands, respectively. Figure 4 shows the behaviour of the SAR time series generated using these invariant samples. This figure describes the monthly mean behaviour of the backscattering, from September 2017 to August 2018, regarding the classes agriculture, forest and grasslands at different channels. In general, these time series

Figure 2. Methodology
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2020, 2020 XXIV ISPRS Congress (2020 edition) present similar patterns. All time series, except for the RGI and CR channels, present a slight increase of the backscattering until December 2017. After that period, they present a decreasing pattern. We assume that this behaviour might be a consequence of changes in precipitation during the referred periods. On the other hand, the RGI and CR time series seem to present a reflected pattern. These time series present opposite behaviours for the three classes over time. The RGI attribute represents a normalized polarization index that theoretically reflects the depolarization thru polarizations and emphasizes the difference between bare ground and forested areas (Hird et al., 2017), which might justify the opposite backscattering pattern observed in this period, in comparison to the other channels. The class separability depends on the temporal class distribution. The ideal channel to use for classification purposes is the one that maximizes the class separability. The CR and RGI time series present high confusion between the backscattering patterns regarding the forest and grasslands classes. On the other hand, these time series show high separability between the class agriculture and the remaining classes, except on April.
Although the VV time series presents temporal backscattering patterns similar to those observed in the VH and NL channels, the VV channel does not present good class separability between January 2018 and March 2018. In this period, for the VV channel, we can observe a confusion between the backscattering patterns regarding areas of agriculture and grasslands. This temporal backscattering similarity between these LULC classes might be associated with the phenological stages of the agricultural cultures observed in the study site. However, extra investigations are required to support this hypothesis.
In general, the temporal backscattering patterns observed in the VH and NL time series indicate higher separability between grasslands and the remaining LULC classes than the other channels considered in this study. Fig. 5 presents the overall accuracies and kappa indices obtained for the classifications performed using different time series. Note that the classification results obtained using the combination of the VV and VH temporal features (described as VV+VH), and the NL features provided similar results. However, the classification performed using VV+VH features provided slightly better results than the other classifications. We performed pairwise comparisons between the kappa indices and the variances of the classification products using the Z-test, in order to verify if there was a significative difference between the results. The Z-test showed that the products obtained using the VV+VH and NL channels are compatible with the hypothesis that there is no significant difference between these products, considering a 95% confidence interval. This test also showed that the products obtained using the NL, VH and VV channels are statistically similar to each other, but slightly different from the results regarding the VV+VH classification, with exception to the NL classification. Moreover, the RGI and CR channels obtained the worst results. The results obtained by these channels were statistically similar, and significantly different to the other classification results. For that reason, in this section, we focus on the results obtained using the VV+VH time series. Figure 6 illustrates the relative importance of each temporal metric estimated during the classification of the VV+VH time series. This figure shows that the most important features for this classification were the coefficient of variation, at VV and VH polarizations, with a relative importance of 20.27% and 19.09%, respectively. The standard deviation was the third most important feature, at the VV polarization, with a relative importance of 7.55%. However, the coefficient of variation represents the ratio of the standard deviation and the mean. Therefore, the temporal metrics related to the standard deviation and the mean might be redundant for the classification process. In this sense, a reduction of the number of temporal metrics may contribute to the classification, with possible improvements to the algorithm's performance and even to the classification accuracy. Figure 3. False-colour composite image that uses temporal metrics to emphasize the difference between areas of agriculture, forest and grasslands, along with the invariant samples Figure 6. Importance of each temporal metric to the RF classification using the combined VV and VH features Figure 7 presents the result of the RF classification of the VV+VH time series. In this classification, we considered all the temporal metrics aforementioned. In addition,  The class Grasslands presented the highest omission error (7%), along with a commission error of 7% as well. In this study, we did not separate the different growing stages of pasture areas into different classes, which might contribute to the obtainment of higher classification errors to the class Grasslands. In addition, grasslands present spectral patterns that might be similar to those observed in areas of agriculture and forest. It might be possible to mitigate these errors with the separation of these patterns into a larger number of classes. We compared our results with the map provide by the MapBiomas project. Figure 8.a illustrates the MapBiomas 2018 reference map, while Figure 8.b represents a difference map between our classification results ( Figure 7) and this reference map. As it shows, the major confusions are related to the borders of the ground features and to small details provided by optical mapping.  b) The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII- B3-2020, 2020XXIV ISPRS Congress (2020 To analyse the quality of the results obtained in this study, we also compared the aforementioned results with results obtained in similar works. Taravat et al. (2019) focused on the development of an automatic grassland cutting status detection based on the analysis of spatio-temporal Sentinel-1 data, using artificial neural networks. The authors obtained classification products with overall accuracy of 85.71%. These results emphasize the potential of machine learning models for grassland cutting status detection with SAR data. Whelen and Siqueira (2018) achieved, over North Dakota, an overall accuracy of 83.5% for grasslands/pasture mapping using Sentinel-1.
Similarly, the authors also reported high confusion between grasslands and agricultural crops, with focus on wheat crops. These results are similar to the results obtained in our study, which also presented confusion between grasslands and agricultural areas.
In addition, in Bazzi et al. (2019) Sentinel-1 TS were used to map paddy rice in France. Although the referred study focused on agricultural mapping, there are similarities to our study, including the use of the RF algorithm to perform the classification of temporal metrics extracted from the TS. However, in this case, the authors used metrics derived from a Gaussian fitting of the VV/VH TS. The adjustment of this methodology for grassland detection might contribute to the reduction of classification confusion observed between areas of grasslands and other vegetation types. The results obtained in the referred study presented an OA of 96.6%. A better understanding of the statistical properties and distribution of the TS data might also contribute to the improvement of the TS classification, as well as testing different speckle filtering techniques.

CONCLUDING REMARKS
This paper explores a Sentinel-1 imagery time series for the classification of grasslands in the northern Brazilian Amazon. We extracted temporal metrics from images acquired over the period of one year, from September 2017 to September 2018, in order to perform a RF classification. The results obtained in this study led us to conclude that the analysis of SAR S1A TS can efficiently assist the detection of grasslands in the northern Brazilian Amazon. We consider that combined VV and VH polarized data are necessary to describe the studied phenomena. The use of TS of indices that aggregate both VV and VH data (NL, CR and RGI) did not result in a significative improvement of the classification, except for the NL index. The z-tests showed that there are no significative differences between the results obtained using the VH+VV data and NL channel.
Future works could also evaluate the potential of using different attributes to analyse each part of the period covered by the time series in order to perform its classification. By doing so, we could create a time series composed of the attributes that best discriminate the analysed classes in each season of the analysed period and also think in the results operationally for mapping the pastures in an automated way.
The results obtained in this study might be complemented with the use of ground data in order to provide more accurate analysis of the classification results, as well as the integration of optical data. The use of ground samples regarding pasture areas with different phenological stages could also assist the detection of different grassland types. In addition, future works can evaluate the use of other classification methods, such as artificial neural network and deep learning techniques, to improve the classification performance.
In addition, it could be interesting to compare the classification results presented in this paper with results generated without the temporal speckle filtering of the TS. The relative calibration of the imagery TS based on a reference image could also mitigate seasonality effects. The application of segmentation techniques could also contribute to the improvement of the classification, due to the complex pixel-based variability that exists because of the speckle effect. Although this study only used attributes derived from the backscatter coefficient, it should be interesting, for future works, the development of methodologies that includes attributes derived from the phase information, making the use of polarimetric and interferometric attributes possible on the discrimination of the LULC classes.