Next Article in Journal
Global Trends in Evapotranspiration Dominated by Increases across Large Cropland Regions
Next Article in Special Issue
A Multi Sensor Approach to Forest Type Mapping for Advancing Monitoring of Sustainable Development Goals (SDG) in Myanmar
Previous Article in Journal
Satellite Observations for Detecting and Forecasting Sea-Ice Conditions: A Summary of Advances Made in the SPICES Project by the EU’s Horizon 2020 Programme
Previous Article in Special Issue
Forest and Land Fires Are Mainly Associated with Deforestation in Riau Province, Indonesia
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Combining Radar and Optical Imagery to Map Oil Palm Plantations in Sumatra, Indonesia, Using the Google Earth Engine

1
Institute of Landscape and Plant Ecology, University of Hohenheim, 70599 Stuttgart, Germany
2
Ecology & Evolutionary Biology, The University of Tennessee, Knoxville, TN 37996, USA
3
National Institute for Mathematical & Biological Synthesis, Knoxville, TN 37996, USA
4
Asian School of the Environment & Earth Observatory Singapore, Nanyang Technological University of Singapore, Singapore 637551, Singapore
*
Author to whom correspondence should be addressed.
Submission received: 2 March 2020 / Revised: 27 March 2020 / Accepted: 8 April 2020 / Published: 10 April 2020

Abstract

:
Monitoring the expansion of commodity crops in the tropics is crucial to safeguard forests for biodiversity and ecosystem services. Oil palm (Elaeis guineensis) is one such crop that is a major driver of deforestation in Southeast Asia. We evaluated the use of a semi-automated approach with random forest as a classifier and combined optical and radar datasets to classify oil palm land-cover in 2015 in Sumatra, Indonesia, using Google Earth Engine. We compared our map with two existing remotely-sensed oil palm land-cover products that utilized visual and semi-automated approaches for the same year. We evaluated the accuracy of oil palm land-cover classification from optical (Landsat), radar (synthetic aperture radar (SAR)), and combined optical and radar satellite imagery (Combined). Combining Landsat and SAR data resulted in the highest overall classification accuracy (84%) and highest producer’s and user’s accuracy for oil palm classification (84% and 90%, respectively). The amount of oil palm land-cover in our Combined map was closer to official government statistics than the two existing land-cover products that used visual interpretation techniques. Our analysis of the extents of disagreement in oil palm land-cover indicated that our map had comparable accuracy to one of them and higher accuracy than the other. Our results demonstrate that a combination of optical and radar data outperforms the use of optical-only or radar-only datasets for oil palm classification and that our technique of preprocessing and classifying combined optical and radar data in the Google Earth Engine can be applied to accurately monitor oil-palm land-cover in Southeast Asia.

Graphical Abstract

1. Introduction

Tropical forests are important ecosystems for biodiversity and carbon storage and they provide numerous material and non-material benefits to people [1]. However, these ecosystems are threatened as a result of fires and conversion for other land-uses, chiefly for the purpose of agriculture [1]. Recent global assessments on the drivers of forest degradation and loss showed that commodity-driven deforestation, which is the loss of forest cover for commodity crops, such as palm oil and soy, was the major driver for deforestation in the tropics [2,3]. This is especially pronounced in tropical Southeast Asia, where one of the major drivers for this land-use change is the expansion of the oil palm (Elaeis guineensis), an economically important vegetable oil crop for the purpose of food, oleochemical industries, and biofuel [4].
Indonesia is the largest producer of palm oil in the world, with 12 million ha of oil palm cultivated and 34 Mt exported in 80 countries, making up 48% of global palm oil export in 2017 [5]. The oil palm industry employs more than 2 million people in Indonesia and contributes to 1.5% of national gross domestic production [6]. However, rapid expansion of oil palm plantations has also led to large-scale deforestation, greenhouse gas emissions, and biodiversity loss, as well as land grabs and social conflicts in Indonesia [7,8,9], leading to many concerns about the negative effects of palm oil on people and the environment [4,10,11,12]. Hence, it is crucial to map and monitor oil palm plantations for evaluating and managing their social and environmental effects [13,14,15].
Multiple techniques have been developed to map oil palm land-cover, using satellite remote sensing data across the tropics [16,17,18,19,20]. Optical methods for mapping oil palm land-cover rely on information derived from the phenology or image characteristics of the oil palm plantation [17]. Phenology-based methods, such as the use of temporal changes in vegetation greenness by Gutierrez-Velez et al. (2011) to detect the area deforested by large-scale oil palm expansion in the Peruvian Amazon. Image-based methods utilize spectral signatures, as well as textural information, to differentiate oil palm trees from their surroundings [21,22]. Oil palm plantations can be manually digitized from satellite images, based on the unique textural information of oil palm plantations (e.g., long rectangular blocks for industrial plantations, geometric shape of oil palm canopy, presence of roads), along with expert knowledge on the land system, or it could be automated based on spectral image analysis, which classifies pixels based on their spectral class thresholds [15]. However, it is challenging to separate oil palm plantations from other spectrally similar land-cover types (e.g., forests, rubber trees) [23], and the frequent presence of cloud cover in the tropics hinders image analysis [17].
To resolve challenges of mapping oil palm plantations in the humid tropics, radar data, which is all-weather and all-time capable, have been tested for mapping oil palm. Microwave energy scattered by vegetation depends on the size, and density, as well as the orientation and dielectric properties of elements that are comparable to the size of the radar wavelength, while optical energy reflected by vegetation depends on the leaf structure, pigmentation, and moisture [24]. Hence, radar data provide more information on the structural properties of the land, such as harvesting trails and roads, which are a prominent feature of oil palm plantations. Radar data are also able to distinguish closed canopy oil palm stands, with their large fronds supported by branchless trunks from other woody vegetation and tree cover [18]. Past studies have successfully used radar data to detect new fronts of deforestation [25], as well as tree-crop plantations in particular oil palm [17,18,26,27]. As a result, a combination of optical and radar data has been used for mapping oil palm across different geographical contexts and using different data combinations [28,29,30,31,32]. In a review of the performance of optical and radar data fusion for land-use mapping and monitoring, 28 out of 32 studies obtained higher accuracies when fusing the two datasets, compared to using only one of them [24]. More specifically for palm oil land cover classification, fusing optical and radar datasets has also resulted in higher accuracies than classification using only one dataset [33].
While the use of optical and radar imagery shows high accuracy for mapping plantation agriculture, automated mapping of oil palm plantations can be problematic and has been shown to underestimate extents of oil palm plantations, as compared to visual classification [15]. In Borneo, automated classification of oil palm was limited to closed canopy oil palm stands, while visual classification was able to detect all areas used for oil palm agriculture, regardless of the phase or condition of the oil palm plantation [15]. This difference is related to automated approaches that are only able to interpret physical characteristics of land-cover, while visual approaches integrate other land-use characteristics, using contextual information and interpreter knowledge of the area into the mapping of oil palm plantations [15].
We mapped oil palm plantations for the year 2015 to test the use of combined radar and optical data to map oil palm plantations in three provinces of Sumatra – Riau, Jambi, and South Sumatra. We compared the accuracy of our combined map with maps generated using only radar or optical data. We also compared the oil palm area of our highest accuracy map with two other oil palm land-cover products produced in 2015: Miettinen et al. (2016), which used Moderate Resolution Imaging Spectroradiometer (MODIS) imagery (250m) to conduct a semi-automated unsupervised land-cover classification with Sentinel-1 data to detect oil palm plantations, and Austin et al. (2017), which manually digitized the extents of oil palm planted using Landsat imagery (30 m) and official plantation boundaries provided by the Indonesian government. Through our study, we aim to address the following research questions: (i) does the use of radar and optical imagery for oil palm classification outperform radar-only and optical-only classification?; (ii) how does our 30 m semi-automated map of oil palm plantations compare with two maps a 250 m semi-automated map and a 30 m digitized map of oil palm plantations for the same year? The study proposes improved methods to map oil palm, using a combination of radar and optical imagery on the Google Earth Engine and a systematic comparison of previous oil palm maps to identify similarities to drive advancement of new methods.

2. Materials and Methods

2.1. Study Site

Our study site were the provinces of Riau (0.5333°N, 101.4500°E), Jambi (1.5833°S, 103.6166°E), and South Sumatra (2.9789°S, 104.7584°E), with a total area of 237,000 km2 or 23,700,000 ha. These three provinces lie east of the Barisan mountain range, spanning an elevation gradient from 0 to 3736 m. Approximately 30% of the land is underlain by peatland close to the eastern coast, with peat depth up to 8 m [34]. These three provinces are also important oil palm producing provinces, with a total growing extent of ~4 million ha, representing 36% of Indonesia’s total oil palm cultivated extent in 2015 [5]. Riau has the largest extent of oil palm cultivated among the three provinces (2,400,876 ha), of which more than half are owned by smallholders (57% or 1,354,503 ha) and the rest (43% or 1,046,373 ha) are owned by industrial government or private estates. South Sumatra has the second largest extent of cultivated oil palm (952,083 ha), of which 45% (431,104 ha) are owned by smallholders and 55% (520,979 ha) are owned by industrial government or private estates. Jambi has the smallest extent of oil palm cultivated (714,399 ha), but the highest proportion of smallholder owned oil palm (63% or 450,075 ha), with the remaining 37% (264,324 ha) owned by industrial government or private estates.

2.2. Satellite Imagery

We chose optical and radar satellite data for the year 2015 to map oil palm plantations in Riau, Jambi, and South Sumatra. Mapping oil palm land-cover for the year 2015 would allow us to compare our extent of oil palm planted area with maps produced for the same year by Miettinen et al. (2016) and Austin et al. (2017), hereafter labelled as ‘Miettinen2016’ and ‘Austin2017’, respectively. Landsat 8 surface reflectance products were imported from the Google Earth Engine data catalogue for our three provinces in Sumatra. Our Landsat imagery consisted of 12 bands, including ultra-blue, blue, green, red, near infrared, and shortwave infrared 1 and 2. Advanced land observing satellite-2 phased array L-band synthetic aperture radar-2 (ALOS-2/PALSAR-2) imagery was imported from the Japan Aerospace Exploration Agency (JAXA) in Google Earth Engine’s data catalogue. Data layers under ALOS-2/PALSAR-2 include horizontal transmit-vertical receive (HV) polarization and horizontal transmit-horizontal receive (HH) polarizations.

2.3. Workflow for Image Processing, Classification, and Analysis

Our workflow includes the following steps: image pre-processing, image classification, accuracy assessment, and comparison with other oil palm maps (Figure 1). Image pre-processing for Landsat imagery removed cloud cover from the surface reflectance Landsat product to create the best available pixel image composite, while that of PALSAR imagery included speckle filtering and conversion of normalized radar cross-sections, as well as the calculation of texture measures and indices. Image classification involved the creation of image stacks using Landsat and PALSAR layers, drawing regions of interest and classifying the image using a machine learning classifier. Three stacks of data layers were created for image classification: PALSAR-only layers (‘synthetic aperture radar (SAR) layer’), Landsat-only layers (‘Landsat layer’), and a combination of Landsat and SAR layers (‘Combined layer’). Accuracy assessments were conducted for SAR, Landsat, and Combined layers to evaluate which map produced the best oil palm extent map. We compared our best oil palm land-cover map based on overall classification accuracy with Miettinen2016 and Austin2017 (1) by examining areas of disagreement against high-resolution Google Earth Pro imagery and (2) by calculating pairwise Jaccard similarity for oil palm cover between the best map produced from our image classification and each of the two existing maps.

2.4. Image Pre-Processing for Landsat Imagery

To remove cloud cover in Landsat images, Landsat image composites were based on a time series, where the best available observation from 40 Landsat 8 images of the path 127 and rows 59 and 60 from the year 2015 were selected through pixel-based image compositing (Figure 2). The image compositing the Google Earth Engine script was adapted from De Alban et al. (2018) and applied to the Landsat 8 surface reflectance product. The image compositing script selected the best available pixel in Landsat 8 images. We applied the median compositing method, which takes the median value of each band and therefore is very resistant to outliers. We used a hybrid cloud shadow masking method, including a function of mask for cloud detection, a cloud score function for cloud masking, and a temporal dark outlier mask (TDOM) for cloud shadow. The threshold for masking clouds was <5% and the shadows z-score = −1. We used the tier 1 surface reflectance Landsat product and Landsat 8 OLI/TIRS sensors.
With the processed Landsat image composite, we calculated five vegetation indices, including the normalized difference vegetation index (NDVI) [35,36] the enhanced vegetation index (EVI) [37,38], the soil adjusted total vegetation index (SATVI) [39], the normalized difference tillage index (NDTI) [40], and the normalized difference water index (NDWI) [41,42], also called the modified NDVI or land surface water index (LSWI; Table 1). The NDVI and LSWI can separate land-cover types, such as forests and croplands [43]. The NDTI, SATVI, and NDWI have been shown to discriminate plantations from forests [32]. Other studies integrating radar and optical satellite imagery have shown that the NDVI, NDWI, and EVI improve the detection of deciduous rubber plantations and forests [32,44,45,46].

2.5. Image Pre-Processing for SAR Imagery

Our 2015 SAR mosaic was obtained from the seamless yearly mosaic available on the Google Earth Engine database and clipped around Sumatra, our area of interest. We applied a 3 × 3 Refined Lee filter on each channel, HH and HV, to smoothen the image, reduce the speckle effects in the raw SAR imagery, and conserve the information of the edges and linear features between two land-cover types [47]. The digital numbers of the HH and HV images were converted into sigma-naught values, following Equation (1).
σ 0 = 10 l o g 10   ( D N 2 ) + CF
σ 0 is the radar backscatter per unit area; DN is the digital number, and CF is the calibration factor with a given value of –83.0 dB for the SAR mosaic [47]. Speckle filtering and conversion to normalized radar cross-section was implemented in the Google Earth Engine based on a script adapted from the SNAP 3.0 toolbox of the European Space Agency (Figure 2).
To complement the information obtained through the SAR channels, HH and HV, we calculated texture measures, which are the relationship between neighbouring pixels and their grey tones in a defined window. The second-order texture measures are the image statistical distribution of grey tone and describe how often a grey tone appears in a specific spatial relationship to another grey tone. To extract textural information, it is necessary to calculate grey-level co-occurrence matrices (GLCMs). GLCM statistic values allow us to understand the organization of pixels and show variation among pixels in different parts of a picture [48]. They outline the angular relationships and distance between adjacent pixels [48]. Texture measures add information to discriminate land-covers and enhance classification accuracy of broad land-cover types [49]. They can predict vegetation properties, such as height, structure, and canopy cover [48].
The GLCM texture measures we calculated included angular second moment (ASM), which represents the orderliness of a picture, the average (AVG), which is the frequency of occurrence of one grey tone in combination with another certain grey tone, the contrast (CON), which is the difference between pixels, the correlation (COR), which is the predictability and linear relationship between two neighbouring pixels, entropy (ENT), which is the degree of disorder among pixels, the inverse difference moment (IDM), which inversely correlated with the contrast, and the variance (VAR), which is the dispersion around the average of pixel values within a defined window. Texture measures were calculated from the average of the directional bands within a 3x3 window using the glcmTexture function available in the Google Earth Engine [50].
To add layers of information and enhance classification, we calculated six indices (Table 2) from the raw synthetic aperture radar backscattering coefficient, including average (AVE), difference (DIF), and HH and HV channels’ simple ratios (RAT1 and RAT2). These four indices yielded good accuracies when differentiating forest to non-forest and coconut palm [51], oil palm, and rubber [33]. In addition, they improve mapping of broad land-cover types, such as cropland, forest, built-up areas, and water [26,45]. The DIF showed a good result to discriminate plantations, such as oil palm and rubber [18]. The two other indices calculated were the normalized difference index (NDI), which is useful to detect new deforestation events [25], and the NL index, which can enhance overall accuracies for land-cover classification using L-band SAR data [52].

2.6. Creating Image Stacks

The final Landsat image stack contained 12 layers, comprising 5 indices (NDVI, EVI, SATVI, NDTI, and LSWI) and 7 bands (blue band (BLUE), green band (GREEN), near-infrared band (NIR), red band (RED), shortwave infrared bands 1 and 2 (SWIR1 and SWIR2), TIR). The final SAR image stack contained 24 layers (i.e., HH and HV with their corresponding textures, AVE, DIF, NDI, NLI, RAT1, and RAT2). The final Combined stack contained 36 layers, including indices, textures, and bands from the two previous stacks. Despite the importance of image features to increase class separability and classification accuracy, increasing the number of features does not necessarily translate in higher accuracies [53]. Each feature was chosen according to its potential to improve oil palm classification, as cited in the literature.

2.7. Regions of Interest and Sampling Design

We defined six land-cover classes, which included forest (FOR), forest shrub mosaic (FSM), oil palm (OPM), bare ground-crop (BGC), built-up area (URB), and water (WTR) (Supplementary Table S1). We delineated regions of interest (ROIs) for these land-cover classes by visual interpretation of high-resolution images in Google Earth Pro from 01 January 2015 to 31 December 2015. We used features such as location, size, shape, tone, colour, shadow, texture, and pattern to identify and differentiate our land-cover classes in Google Earth Pro [54] (Supplementary Table S1). To ensure representability of our sample ROIs, we obtained between one and four ROIs for each land-cover class in each regency for each of our provinces of interest; Riau, Jambi, and South Sumatra. This sampling scheme ensured a spread of ROIs across our different satellite image tiles. The final sample of 788 polygons was imported as a feature collection in the Google Earth Engine.
We used the randomColumn and sampleRegions functions in the Google Earth Engine to separate the ROIs into testing and training polygons. The randomColumn function gives a random number between 0 and 1 to each polygon, which would be appointed as training data when they have a value of ≤0.70 and as testing data when they have a value of >0.70. The same selection process was repeated at the pixel level. Using the same Google Earth Engine functions, every pixel within each training and testing polygon was assigned a random number between 0 and 1. Pixels from testing polygons with an assigned value of ≤0.70 were segregated as testing pixels, while pixels from training polygons with an assigned value of ≤0.70 were appointed as training pixels. In total, the sample for the SAR classification had 191,663 training pixels and 108,816 testing pixels, the Landsat sample had 133,064 training pixels and 75,502 testing pixels, and the Combined sample had 133,039 training pixels and 75,501 testing pixels. This two-level random sampling design prevented training and testing data to be spatially dependent [55].

2.8. Image Classification

We used a non-parametric random forest (RF) machine learning classifier for supervised land-cover classification. RF is an ensemble of decision trees that, based on the data they receive, vote for the most popular class for an input pixel [56]. In general, non-parametric machine-learning classifiers, such as RF or support vector machine (SVM), yield higher accuracies than parametric classifiers. RF can improve accuracy of land-cover classification [57] and detect plantations like oil palm [54]. We used the RandomForest function included in the Google Earth Engine and trained the machine learning algorithm with the training data obtained from each of our three image stacks: Landsat, SAR, and Combined. The following parameters were used as inputs for the algorithm: 100 for the number of decision trees created per class, 10 for the minimal length of a terminal node, and 0.50 as the input fraction to bag per tree. Following classification, we applied a majority filter to homogenise the classified image and remove small clusters of pixels. Isolated pixels were replaced by another pixel value based on the majority of their contiguous neighbouring pixels.

2.9. Accuracy Assessments

We evaluated our classification results using standard accuracy assessment metrics as the error matrix, overall accuracy, user’s, and producer’s accuracies [58]. Since we were more interested in oil palm mapping, we focused on the producer’s and user’s accuracies for oil palm land-cover class. We identified the best map as the one with the highest overall, producer’s, and user’s accuracies. To calculate the metrics above, we used Google Earth Engine functions errorMatrix(), accuracy(), consumersAccuracy(), producersAccuracy(), and kappa().

2.10. Comparing Classified Oil Palm Extent with Other Oil Palm Maps for 2015

We extracted the oil palm mapped under our best classified map and compared this with two other oil palm land-cover products produced in 2015 from Miettinen2016 and Austin2017 using difference images. Our difference images allowed for the visualization of areas of agreement, where both maps detected oil palm or both maps detected no oil palm, as well as areas of disagreement, where one map detected oil palm but not the other [59]. We first used the function r.cross from GRASS GIS to create the two difference images, comparing our classified map with Miettinen2016 and Austin2017. The raster layers were polygonised to calculate the areas of agreement and disagreement in each map comparison. These values of agreement and disagreement were used to calculate the Jaccard similarity coefficient (J) using the following equations:
J ( X , Y ) = | X Y | / | X Y |
Under Equation (2), the Jaccard coefficient (J) is calculated by dividing the areas shared by maps X and Y by the areas belonging to both or either maps. The Jaccard coefficient values vary between 0 and 1. A value of 0 means that the two compared maps have no areas in common, while a value of 1 means that the two maps are identical and share the same area.
To investigate how land-cover classes were classified differently under our extents of disagreement between our classified oil palm extent and Miettinen2016, as well as our classified oil palm extent and Austin2017, we extracted the areas of disagreement in each of our difference images. These areas of disagreement refer to pixels identified as oil palm in one map extent but not the other (i.e., when our map identified oil palm but Miettinen2016 or Austin2017 did not and where Miettinen2016 or Austin2017 identified oil palm where our map did not). For each extent of disagreement, we randomly selected 100 points.
We visually inspected each of these points in Google Earth Pro, using high-resolution images from 2015, and recorded the land-cover for each point. We counted the number of points for each land-cover class and grouped them according to whether the point was identified as oil palm under our map or either Miettinen2016 or Austin2017. Since the area of oil palm identified under Miettinen2016, Austin2017, and our maps differed, the random sample of 100 points from our disagreement extents reflected this difference in area. Hence, we weighted our count data of each land-cover class with the total number of random points drawn from each map. For example, 56 pixels were correctly identified as oil palm in the Combined map and 15 in Miettinen2016; however, after weighting these count data, both maps had 71% of their pixels correctly identified as oil palm. We present the proportion of points that were accurately identified as oil palm in our Google Earth Pro image and the proportion of points which fell under other land-cover classes, such as forest shrub mosaic or bare ground cover land-classes.

3. Results

3.1. Comparing Landsat, SAR, and Combined Maps

Our Combined (Landsat + SAR) map had the highest overall accuracy score of 84%, followed by our Landsat map (78%) and SAR map (74%) (Table 3, Figure 3). The Combined map also had the highest producer’s and user’s accuracies for oil palm land-cover class (84% and 90%, respectively), as compared to Landsat (70% and 62%) and SAR (81% and 66%) maps (Table 3). The error matrix for our Landsat classified image showed a high number of oil palm pixels misclassified as forest and forest shrub mosaic land-cover classes (Figure 4). Under our SAR classified image, oil palm pixels were confused with mostly bare ground-crop, followed by forest shrub mosaic and forest land-cover classes (Figure 4). The error matrix for our Combined classified image showed that oil palm pixels were confounded with bare ground-crop and forest shrub mosaic (Figure 4). Overall, combining SAR and Landsat imageries outperformed single imagery (SAR-only and Landsat-only) classification in mapping oil-palm cover due to the superior producer ‘s and user’s accuracy of the former technique (Table 3).

3.2. Comparing Oil Palm Extents with Miettinen2016 and Austin2017

The extent of oil palm mapped in our Combined classified image for Riau, Jambi, and South Sumatra was 4,418,176 ha, which was 8.6% higher than the oil palm extent from the Indonesian Ministry of Tree Crop Statistics (4,067,358 ha) for the year 2015 [5]. Our mapped extent of oil palm was 47% higher than the extent of oil palm mapped by Miettinen2016 (2,997,943 ha) and 22% higher than the extent of oil palm mapped by Austin2017 (3,608,815 ha).
The Jaccard coefficient between our Combined map and the Miettinen2016 map was 0.55, which was higher than the Jaccard coefficient between our Combined map and the Austin2017 map (0.40).
The difference images in Figure 5 show extents of agreement and disagreement for the Combined-Miettinen2016 and Combined-Austin2017 maps. These difference images show the locations of agreements, where oil palm was detected in both maps, and locations of disagreements, where oil palm was detected in one map but not the other. In the Combined-Miettinen2016 comparison, the total extent of disagreement was 2,155,972 ha, of which 83% (1,788,479 ha) stemmed from our Combined map, identifying oil palm in areas where Miettienen2016 did not, and 17% (367,496 ha) were from areas of disagreement where Miettinen2016 identified oil palm and where our Combined map did not. In the Combined-Austin2017 comparison, the total extent of disagreement was 3,444,846 ha, of which 62% (2,126,000 ha) came from our Combined map, identifying oil palm in areas where Austin2017 did not, and 38% (1,318,846 ha) were from areas of disagreement where Austin2017 identified as oil palm and where our Combined map did not. These extents were largely concentrated in South Sumatra. Approximately 62% (2,125,999 ha) of the extent of disagreement was a result of our Combined map detecting oil palm where Austin2017 did not, and these areas were found mostly on the eastern coasts of Riau and Jambi.
Of the 100 random sample points for our Combined-Miettinen2016 comparison, 79 points were drawn from the Combined map and 21 points were drawn from the Miettinen2016 map. Whereas, for our Combined-Austin2017 map, 58 points were drawn from the Combined map, while 42 points were drawn from the Austin2017 map. We calculated the weighted percentage of points for each land-cover class using the total number of points identified under each map above (Figure 6).
The majority of the points identified as oil palm in Combined (70.8%) and Miettinen2016 (71.4%) were verified as oil palm under Google Earth Pro (Figure 6). About 15.2% of the oil palm points under the Combined map were forest shrub mosaic, and 13.9% of these points were identified as bare ground crop with Google Earth Pro. For the Miettinen2016 map, 19% of the oil palm points were forest shrub mosaic, 4.8% were bare ground cover, and another 4.8% identified as urban extents under Google Earth Pro (Figure 6). Under the Combined-Austin2017 comparison, the majority of the points identified as oil palm in Combined (65.5%) were verified as oil palm under Google Earth Pro, while 28.5% of oil palm points under Austin2017 were verified as oil palm under Google Earth Pro (Figure 6). About 38.1% of the oil palm points under Austin2017 were forest shrub mosaic, 28.6% were bare ground cover, and 4.8% were forest under Google Earth Pro. Under the Combined map, 22.4% of the oil palm points were bare ground cover, while 12.1% were forest shrub mosaic under Google Earth Pro (Figure 6).

4. Discussion

Our results showed that combining optical and radar datasets improved the accuracy of oil palm land-cover classification. This corroborates with results from other recent studies that use both optical and radar data to map oil palm plantations in Myanmar [32,33] and Indonesia [30]. Moreover, we developed a code to pre-process radar and optical datasets, as well as classify land-cover entirely in the Google Earth Engine. The improved accuracy for detecting oil palm plantations with the use of optical and radar data highlights the importance of multisensory and data fusion techniques, which integrate multiple data features to improve map detail and accuracy [60]. This is especially important in the tropics, where oil palm plantations expand most and where optical-based mapping is limited due to high extents of cloud cover. In our study, the use of only Landsat data underestimated the extent of oil palm plantations due to confusion with spectrally similar land-cover classes, such as forest and forest shrub mosaic, while the use of only SAR data overestimated the extent of oil palm plantations, largely as a result of confusion between bare ground cover and oil palm (Figure 4). In general, classification using SAR data overlooks new oil palm plantations, because they have similar backscatter as herbaceous vegetation and bare land [17].
Our combined Landsat and SAR approach produced an oil palm extent of 4,418,176 ha, which was closest to the oil palm area reported by the Indonesian Ministry of Tree Crop Statistics (4,067,358 ha), as compared to Miettinen2016 (2,997,943 ha) and Austin2017 (3,608,815 ha) for the provinces of Riau, Jambi, and South Sumatra. However, disagreement maps between our combined approach and Miettinen2016 and Austin2017 showed large extents of non-overlap in the identification of oil palm land-cover, with our Combined map showing a higher level of agreement with the Miettinen2016 map as compared to the Austin2017 map. These differences are mainly the result of different methods used to define and classify oil palm land-cover in these studies.
Our semi-automated supervised 30 m Combined map was at a higher resolution compared to Miettinen2016, which used Sentinel-1 C-band radar data (20 m) that was resampled to 250 m to match MODIS images [30]. This difference in spatial resolution and Miettinen2016’s definition of palm plantations as contiguous closed canopy plantations larger than 1 km2 could be one reason for the difference in oil palm area (Combined map had 47% larger extent), as well as the differences in misclassification, where a higher proportion of oil palm points were forest land-cover under Miettinen2016, as compared to a higher proportion of oil palm points as bare ground cover under the Combined map (Figure 6). While some pixels from Miettinen2016 could be coconut and sago palm plantations, the authors visually inspected and removed any potential misclassified areas to avoid classification errors. This was not carried out in our study. Given that coconut and sago plantations are prominent along the coasts of Sumatra (Lee J.S.H. pers. obs.), the large extents of disagreement between our Combined map and Miettinen2016 on the east coast of Sumatra (Figure 5A) are likely coconut and sago plantations misclassified as oil palm in our Combined map. Our results support previous studies that found higher accuracies using supervised methods [61,62,63].
Austin2017 manually digitized oil palm by visually interpreting Landsat imagery, using rectilinear patterns and context indicators, such as road networks, to identify oil palm plantations [64]. Only large-scale oil palm plantations that were organized in a grid pattern and associated with physical structures, such as roads and mills, were mapped. In addition, recently cleared extents proximal to existing plantations were included. The focus of Austin2017 on large-scale oil palm plantations and exclusion of small-scale plantations could be the reason why the Combined map detected a larger extent of oil palm area (22% higher than Austin2017). The focus on land-use by Austin2017’s visual interpretation approach is an alternative to automated detection methods, such as our Combined map, which limits the interpretation of physical data to that of land-cover [15]. Hence, the differences in misclassifications of land-cover under Google Earth Pro, where Austin2017 identified a higher number of forest shrub mosaic and bare ground cover as oil palm, may not be so much a misclassification of oil palm land-cover under Austin2017’s approach, but rather may reflect the different development phases of the plantation or the condition of the crop (e.g., immature, diseased) due to Austin2017’s use of auxiliary information and knowledge of the area. This difference was also demonstrated using automated and visual interpretation approaches to oil palm mapping in Borneo by Miettinen et al. (2017) and Gaveau et al. (2016), respectively [15].

5. Conclusions

Our study shows the usefulness of the Google Earth Engine for processing large satellite datasets and producing classification results of land-cover change at regional scales. Specifically, our study highlights the advantages of combining radar and optical data to improve the accuracy of oil palm land-cover classification. The accessibility of optical and radar datasets on the Google Earth Engine, along with publicly available documentation of code, allows for collaboration and cooperation among different users to conduct land-use and land-cover change for environmental monitoring [65]. Our code for pre-processing Landsat and SAR data and conducting land-cover classification is available under the Supplementary Methods and can be used to improve semi-automated techniques to map and monitor oil palm expansion across the tropics [66]. However, our study also highlights limitations of automated mapping techniques for oil palm plantations and seconds the conclusion by Miettinen et al. (2018) about the need to combine use of visual and automated oil palm mapping approaches for monitoring oil palm expansion. Contextual factors, such as the presence of other palm plantations and the phase of development of oil palm plantations, are some examples of information that cannot be inferred using satellite imagery alone and revealed the differences between our combined optical and radar automated approach and mixed visual and automated mapping approaches from Miettinen2016 and Austin2017.

Supplementary Materials

The following are available online at https://0-www-mdpi-com.brum.beds.ac.uk/2072-4292/12/7/1220/s1: Table S1: Description of 5 land-cover classes used in image classification; Table S2: Number of training and testing pixels for each image stack; Table S3: Error matrices for Landsat classified image; Table S4: Error matrices for SAR classified image; Table S5: Error matrices for Combined classified image; Table S6: Identification reference of the 40 Landsat images used in the classification; Table S7: Producer’s (%) and user’s (%) accuracy of Combined, Landsat, and SAR land-cover maps, classified with a Classification and Regression Tree (CART) algorithm. Code for data pre-processing and land-cover classification is available in Github repository. (https://github.com/thuansarzynski/GEE_CombinedLandsatSAR_oilpalm)

Author Contributions

Conceptualization of the research project and methods involved T.S., X.G., L.C., and J.S.H.L. T.S. was primarily responsible for developing the code used in this analysis. Formal analysis, validation of the maps, and data curation was conducted by T.S. Writing and draft preparation was carried out by T.S. and J.S.H.L., while writing, reviewing, and editing was done by X.G. and L.C. Data visualization was carried out by T.S. and J.S.H.L. Supervision of T.S. was from J.S.H.L. Project administration and funding was administered by J.S.H.L. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by Singapore Ministry of Education Academic Research Fund Tier 1 grant number RG146/16.

Acknowledgments

The authors wish to thank the Google Earth Engine community and especially Noel Gorelick and Pablo Martin, who consistently answered our questions and debugged our codes. The authors thank Jose Don De Alban for reviewing the manuscript and providing useful comments. The authors also thank Pern Zee Wei and Nicholas Ang Kang Jing, who participated in data collection.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Díaz, S.; Settele, J.; Brondízio, E.; Ngo, H.T.; Guèze, M.; Agard, J.; Arneth, A.; Balvanera, P.; Butchart, S.; Chan, K.; et al. Intergovernmental Platform on Biodiversity and Ecosystem Services (IPBES). Available online: https://ipbes.net/global-assessment (accessed on 9 April 2020).
  2. Curtis, P.G.; Slay, C.M.; Harris, N.L.; Tyukavina, A.; Hansen, M.C. Classifying drivers of global forest loss. Science 2018, 361, 1108–1111. [Google Scholar] [CrossRef] [PubMed]
  3. Vijay, V.; Reid, C.D.; Finer, M.; Jenkins, C.N.; Pimm, S.L. Deforestation risks posed by oil palm expansion in the Peruvian Amazon. Environ. Res. Lett. 2018, 13, 110–118. [Google Scholar] [CrossRef]
  4. Meijaard, E.; Garcia-Ulloa, J.; Sheil, D.; Wich, S.A.; Carlson, K.M.; Juffe-Bignoli, D.; Brooks, T.M. Oil Palm and Biodiversity. A Situation Analysis by the IUCN Oil Palm Task Force; IUCN: Gland, Switwerland, 2018. [Google Scholar]
  5. Hendaryati Indonesian Oil Palm Statistics. In Indonesian Oil Palm Statistics; Indonesian Ministry of Agriculture: Jakarta, Indonesia, 2017.
  6. Mukherjee, I.; Sovacool, B.K. Palm oil-based biofuels and sustainability in Southeast Asia: A review of Indonesia, Malaysia, and Thailand. Renew. Sustain. Energy Rev. 2014, 37, 1–12. [Google Scholar] [CrossRef]
  7. Nilakrisna, D.; Taqwa, R.; Purnama, D.H.; Susanto, R.H. Social Conflicts Between Oil-Palm Plantation Company and Indigenous People in Jambi Province. J. Trop. Life Sci. 2016, 6, 113–117. [Google Scholar]
  8. Abood, S.A.; Lee, J.S.H.; Burivalova, Z.; Garcia-Ulloa, J.; Koh, L.P. Relative Contributions of the Logging, Fiber, Oil Palm, and Mining Industries to Forest Loss in Indonesia. Conserv. Lett. 2015, 8, 58–67. [Google Scholar] [CrossRef] [Green Version]
  9. Meijaard, E.; Sheil, D. Oil-Palm Plantations in the Context of Biodiversity Conservation; Elsevier Academic Press: Amsterdam, The Netherlands, 2013; Volume 5, ISBN 9780123847195. [Google Scholar]
  10. Rist, L.; Feintrenie, L.; Levang, P. The livelihood impacts of oil palm: Smallholders in Indonesia. Biodivers. Conserv. 2010, 19, 1009–1024. [Google Scholar] [CrossRef]
  11. Forsyth, T. Public concerns about transboundary haze: A comparison of Indonesia, Singapore, and Malaysia. Glob. Environ. Chang. 2014, 25, 76–86. [Google Scholar] [CrossRef] [Green Version]
  12. Austin, K.G.; Kasibhatla, P.S.; Urban, D.L.; Stolle, F.; Vincent, J. Reconciling oil palm expansion and climate change mitigation in Kalimantan, Indonesia. PLoS ONE 2015, 10, e0127963. [Google Scholar] [CrossRef]
  13. Gaveau, D.L.A.; Sheil, D.; Husnayaen; Salim, M.A.; Arjasakusuma, S.; Ancrenaz, M.; Pacheco, P.; Meijaard, E. Rapid conversions and avoided deforestation: Examining four decades of industrial plantation expansion in Borneo. Sci. Rep. 2016, 6, 1–13. [Google Scholar] [CrossRef]
  14. Vijay, V.; Pimm, S.L.; Jenkins, C.N.; Smith, S.J. The Impacts of Oil Palm on Recent Deforestation and Biodiversity Loss. PLoS ONE 2016, 11, e0159668. [Google Scholar] [CrossRef]
  15. Miettinen, J.; Gaveau, D.L.A.; Liew, S.C. Comparison of visual and automated oil palm mapping in Borneo. Int. J. Remote Sens. 2018, 40, 8174–8185. [Google Scholar] [CrossRef]
  16. Gutiérrez-Vélez, V.H.; DeFries, R.; Pinedo-Vásquez, M.; Uriarte, M.; Padoch, C.; Baethgen, W.; Fernandes, K.; Lim, Y. High-yield oil palm expansion spares land at the expense of forests in the Peruvian Amazon. Environ. Res. Lett. 2011, 6, 44029. [Google Scholar] [CrossRef]
  17. Li, L.; Dong, J.; Tenku, S.N.; Xiao, X. Mapping Oil Palm Plantations in Cameroon Using PALSAR. Remote Sens. 2015, 7, 1206–1224. [Google Scholar] [CrossRef] [Green Version]
  18. Miettinen, J.; Liew, S.C. Separability of insular Southeast Asian woody plantation species in the 50 m resolution ALOS PALSAR mosaic product. Remote Sens. Lett. 2011, 7, 58–63. [Google Scholar] [CrossRef]
  19. Shafri, H.Z.M.; Hamdan, N.; Saripan, M.I. Semi-automatic detection and counting of oil palm trees from high spatial resolution airborne imagery. Int. J. Remote Sens. 2011, 32, 2095–2115. [Google Scholar] [CrossRef]
  20. Srestasathiern, P.; Rakwatin, P. Oil palm tree detection with high resolution multi-spectral satellite imagery. Remote Sens. 2014, 6, 9749–9774. [Google Scholar] [CrossRef] [Green Version]
  21. Carlson, K.M.; Curran, L.M.; Asner, G.P.; Pittman, A.M.D.; Trigg, S.N.; Marion Adeney, J. Carbon emissions from forest conversion by Kalimantan oil palm plantations. Nat. Clim. Chang. 2013, 3, 283–287. [Google Scholar] [CrossRef]
  22. Thenkabail, P.S.; Stucky, N.; Griscom, B.W.; Ashton, M.S.; Diels, J.; van der Meer, B.; Enclona, E. Biomass estimations and carbon stock calculations in the oil palm plantations of African derived savannas using IKONOS data. Int. J. Remote Sens. 2004, 25, 5447–5472. [Google Scholar] [CrossRef]
  23. Morel, A.C.; Fisher, J.B.; Malhi, Y. Evaluating the potential to monitor aboveground biomass in forest and oil palm in Sabah, Malaysia, for 2000–2008 with Landsat ETM + and ALOS-PALSAR. Int. J. Remote Sens. 2012, 33, 3614–3639. [Google Scholar] [CrossRef]
  24. Joshi, N.; Baumann, M.; Ehammer, A.; Fensholt, R.; Grogan, K.; Hostert, P.; Rudbeck Jepsen, M.; Kuemmerle, T.; Meyfroidt, P.; Mitchard, E.T.A.; et al. A review of the Application of Optical and Radar Remote Sensing Data Fusion to Land Use Mapping and Monitoring. Remote Sens. 2016, 8, 70. [Google Scholar] [CrossRef] [Green Version]
  25. Filho Almeida, R.; Shimabukuro, Y.E.; Rosenqvist, A.; Sanchez, G.A. Using dual—Polarized ALOS PALSAR data for detecting new fronts of deforestation in the Brazilian Amazônia. Int. J. Remote Sens. 2009, 30, 37–41. [Google Scholar] [CrossRef]
  26. Dong, X.; Quegan, S.; Yumiko, U.; Hu, C.; Zeng, T. Feasibility Study of C- and L-band SAR Time Series Data in Tracking Indonesian Plantation and Natural Forest Cover Changes. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2015, 8, 3692–3699. [Google Scholar] [CrossRef]
  27. Ramdani, F. Recent expansion of oil palm plantation in the most eastern part of Indonesia: Feature extraction with polarimetric SAR. Int. J. Remote Sens. 2019, 40, 7371–7388. [Google Scholar] [CrossRef]
  28. Cheng, Y.; Yu, L.; Xu, Y.; Liu, X.; Lu, H.; Cracknell, A.P. Towards global oil palm plantation mapping using remote-sensing data. Int. J. Remote Sens. 2018, 39, 5891–5906. [Google Scholar] [CrossRef]
  29. Cheng, B.; Furtado, A.; Smyth, H.E.; Henry, R.J. Influence of genotype and environment on coffee quality. Trends Food Sci. Technol. 2016, 57, 20–30. [Google Scholar] [CrossRef] [Green Version]
  30. Miettinen, J.; Shi, C.; Liew, S.C. 2015 Land cover map of Southeast Asia at 250 m spatial resolution. Remote Sens. Lett. 2016, 7, 701–710. [Google Scholar] [CrossRef]
  31. Miettinen, J.; Shi, C.; Liew, S.C. Towards automated 10–30 m resolution land cover mapping in insular South-East Asia. Geocarto Int. 2017, 34, 443–457. [Google Scholar] [CrossRef]
  32. Torbick, N.; Ledoux, L.; Salas, W.; Zhao, M. Regional Mapping of Plantation Extent Using Multisensor Imagery. Remote Sens. 2016, 8, 236. [Google Scholar] [CrossRef] [Green Version]
  33. DeAlban, J.D.; Connette, G.M.; Webb, E.L. Combined Landsat and L-Band SAR Data Improves Land Cover Classification and Change Detection in Dynamic Tropical Landscapes. Remote Sens. 2018, 10, 306. [Google Scholar] [CrossRef] [Green Version]
  34. Wahyunto; Suryadiputra, N. Peatland Distribution in Sumatra and Kalimantan-Explanation of Its Data Sets Including Source of Information, Accuracy, Data Constraints and Gaps; Wetlands International: Indonesia Programme, Bogor, 2008. [Google Scholar]
  35. Rouse, W.; Haas, H.; Deering, W. Monitoring vegetation systems in the great plains with erts. NASA. Goddard Space Flight Center 1974, 1, 309–317. [Google Scholar]
  36. Tucker, C.J. Red and Photographic Infrared Linear Combinations for Monitoring Vegetation. Remote Sens. Environ. 1979, 150, 127–150. [Google Scholar] [CrossRef] [Green Version]
  37. Huete, A.R.; Liu, H.Q.; Batchily, K.; van Leeuwen, W. A Comparison of Vegetation Indices over a Global Set of TM Images for EOS-MODIS. Remote Sens. Environ. 1997, 59, 440–451. [Google Scholar] [CrossRef]
  38. Huete, A.; Didan, K.; Miura, T.; Rodriguez, E.P.; Gao, X.; Ferreira, L.G. Overview of the radiometric and biophysical performance of the MODIS vegetation indices. Remote Sens. Environ. 2002, 83, 195–213. [Google Scholar] [CrossRef]
  39. Marsett, R.C.; Qi, J.; Heilman, P.; Biedenbender, S.H.; Watson, M.C.; Amer, S.; Weltz, M.; Goodrich, D.; Marsett, R. Remote Sensing for Grassland Management in the Arid Southwest. Rangel. Ecol. Manag. 2006, 59, 530–540. [Google Scholar] [CrossRef]
  40. Van Deventer, A.P.; Ward, A.D.; Gowda, P.H.; Lyon, J.G. Using Thematic Mapper Data to Identify Contrasting Soil Plains and Tillage Practices. Photogramm. Eng. Remote Sens. 1997, 63, 87–93. [Google Scholar]
  41. Gao, B.C. A normalized difference water index for remote sensing of vegetation liquid water from space. SPIE 1996, 2480, 225–236. [Google Scholar] [CrossRef]
  42. Jurgens, C. International Journal of Remote The modified normalized difference vegetation index (mNDVI) a new index to determine frost damages in agriculture based on Landsat TM data. Int. J. Remote Sens. 1997, 18, 3583–3594. [Google Scholar] [CrossRef]
  43. Xiao, X.; Boles, S.; Frolking, S.; Salas, W.; Iii, B.M.; Li, C.; He, L.; Zhao, R. Landscape-scale characterization of cropland in China using Vegetation and Landsat TM images. Int. J. Remote Sens. 2002, 23, 3579–3594. [Google Scholar] [CrossRef]
  44. Chen, B.; Li, X.; Xiao, X.; Zhao, B.; Dong, J. Mapping tropical forests and deciduous rubber plantations in Hainan Island, China by integrating PALSAR 25-m and multi-temporal Landsat images. Int. J. Appl. Earth Obs. Geoinf. 2016, 50, 117–130. [Google Scholar] [CrossRef]
  45. Dong, J.; Xiao, X.; Sheldon, S.; Biradar, C.; Xie, G. Mapping tropical forests and rubber plantations in complex landscapes by integrating PALSAR and MODIS imagery. ISPRS J. Photogramm. Remote Sens. 2012, 74, 20–33. [Google Scholar] [CrossRef]
  46. Kou, W.; Xiao, X.; Dong, J.; Gan, S.; Zhai, D.; Zhang, G.; Qin, Y.; Li, L. Mapping Deciduous Rubber Plantation Areas and Stand Ages with PALSAR and Landsat Images. Remote Sens. 2015, 7, 1048–1073. [Google Scholar] [CrossRef] [Green Version]
  47. Haralick, R.M.; Shanmugam, K.; Dinstein, I. Textural Features for Image Classification. IEEE Trans. Syst. Man Cybern. 1973, 3, 610–621. [Google Scholar] [CrossRef] [Green Version]
  48. Gallardo-cruz, J.A.; Lebrija-trejos, E.E.; Marco, A.; Meave, J.A.; Gonza, E.J.; Herna, L.; Gallardo-cruz, R.; Pe, E.A.; Martorell, C. Predicting Tropical Dry Forest Successional Attributes from Space: Is the Key Hidden in Image Texture? PLoS ONE 2012, 7, e30506. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  49. Rakwatin, P.; Longépé, N.; Isoguchi, O. International Journal of Remote Using multiscale texture information from ALOS PALSAR to map tropical. Int. J. Remote Sens. 2012, 33, 7727–7746. [Google Scholar] [CrossRef]
  50. Conners, R.W.; Trivedi, M.M.; Harlow, C.A. Segmentation of a High-Resolution Urban Scene Using Texture Operators *. Comput. Vis. Graph. Image Process. 1984, 25, 273–310. [Google Scholar] [CrossRef]
  51. Estomata, M.T. Forest Cover Classification and Change Detection Analysis Using ALOS PALSAR Mosaic Data to Support the Establishment of a Pilot MRV System for REDD-Plus on Leyte Island; Deutsche Gesellschaft für Internationale Zusammenarbeit (GIZ) GmbH: Quezon City, Philippines, 2014; pp. 1–18. [Google Scholar]
  52. Li, G.; Lu, D.; Moran, E.; Dutra, L.; Batistella, M. A comparative analysis of ALOS PALSAR L-band and RADARSAT-2 C-band data for land-cover classification in a tropical moist region. ISPRS J. Photogramm. Remote Sens. 2012, 70, 26–38. [Google Scholar] [CrossRef] [Green Version]
  53. Giri, C. Remote Sensing of Land Use and Land Cover; CRC Press: Boca Raton, FL, USA, 2012; ISBN 9781420070750. [Google Scholar]
  54. Lee, J.S.H.; Wich, S.; Widayati, A.; Koh, L.P. Detecting industrial oil palm plantations on Landsat images with Google Earth Engine. Remote Sens. Appl. Soc. Environ. 2016, 4, 219–224. [Google Scholar] [CrossRef] [Green Version]
  55. Zhen, Z.; Quackenbush, L.J.; Stehman, S.V.; Zhang, L. Impact of training and validation sample selection on classification accuracy and accuracy assessment when using reference polygons in object-based classification. Int. J. Remote 2013, 34, 6914–6930. [Google Scholar] [CrossRef]
  56. Breiman, L. Random Forests. Mach. Learn. 2001, 45, 5–32. [Google Scholar] [CrossRef] [Green Version]
  57. Rodriguez-galiano, V.F.; Ghimire, B.; Rogan, J.; Chica-olmo, M.; Rigol-sanchez, J.P. An assessment of the effectiveness of a random forest classifier for land-cover classification. ISPRS J. Photogramm. Remote Sens. 2012, 67, 93–104. [Google Scholar] [CrossRef]
  58. Congalton, R.G.; Green, K. Assessing the Accuracy of Remotely Sensed Data; CRC Press: Boca Raton, FL, USA, 2006; ISBN 9781420055122. [Google Scholar]
  59. Congalton, R.G. Accuracy assessment and validation of remotely sensed and other spatial information. Int. J. Wildland Fire 2001, 10, 321–328. [Google Scholar] [CrossRef] [Green Version]
  60. Nomura, K.; Mitchard, E.T.A.; Patenaude, G.; Bastide, J.; Oswald, P.; Nwe, T. Oil palm concessions in southern Myanmar consist mostly of unconverted forest. Sci. Rep. 2019, 9, 11931. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  61. Mohammady, M.; Moradi, H.R.; Zeinivand, H.; Temme, A.J.A.M. A comparison of supervised, unsupervised and synthetic land use classification methods in the north of Iran. Int. J. Environ. Sci. Technol. 2015, 12, 1515–1526. [Google Scholar] [CrossRef] [Green Version]
  62. Krishna Bahadur, K.C. Improving landsat and irs image classification: Evaluation of unsupervised and supervised classification through band ratios and dem in a mountainous landscape in Nepal. Remote Sens. 2009, 1, 1257–1272. [Google Scholar] [CrossRef] [Green Version]
  63. Ismail, M.H. Evaluating supervised and unsupervised techniques for land cover mapping using remote sensing data. Geog. Malays. J. Soc. Space 2009, 5, 1–10. [Google Scholar]
  64. Austin, K.G.; Mosnier, A.; Pirker, J.; McCallum, I.; Fritz, S.; Kasibhatla, P.S. Shifting patterns of oil palm driven deforestation in Indonesia and implications for zero-deforestation commitments. Land Use Policy 2017, 69, 41–48. [Google Scholar] [CrossRef] [Green Version]
  65. Gorelick, N.; Hancher, M.; Dixon, M.; Ilyushchenko, S.; Thau, D.; Moore, R. Google Earth Engine: Planetary-scale geospatial analysis for everyone. Remote Sens. Environ. 2017, 202, 18–27. [Google Scholar] [CrossRef]
  66. Koh, L.P.; Miettinen, J.; Liew, S.C.; Ghazoul, J. Remotely sensed evidence of tropical peatland conversion to oil palm. Proc. Natl. Acad. Sci. USA 2011, 108, 5127–5132. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Overall workflow for mapping oil palm plantations in Riau, Jambi, and South Sumatra, which includes pre-processing of our phased array L-band synthetic aperture radar-2 (PALSAR-2) and Landsat 8 images, image classification using a random forest algorithm, developing accuracy assessments for our oil palm land cover map, and comparing our oil palm land cover map with two other oil palm map sources.
Figure 1. Overall workflow for mapping oil palm plantations in Riau, Jambi, and South Sumatra, which includes pre-processing of our phased array L-band synthetic aperture radar-2 (PALSAR-2) and Landsat 8 images, image classification using a random forest algorithm, developing accuracy assessments for our oil palm land cover map, and comparing our oil palm land cover map with two other oil palm map sources.
Remotesensing 12 01220 g001
Figure 2. Results of image pre-processing for synthetic aperture radar (SAR) and Landsat imagery taken from the regency of Kampar, Riau province, 101°25′00″E 0.50°00′00″N.
Figure 2. Results of image pre-processing for synthetic aperture radar (SAR) and Landsat imagery taken from the regency of Kampar, Riau province, 101°25′00″E 0.50°00′00″N.
Remotesensing 12 01220 g002
Figure 3. Classified image for Riau, Jambi, and South Sumatra using our Combined image stack. Insets focus on the border of Riau and Jambi for (A) Combined, (B) SAR, and (C) Landsat classified images.
Figure 3. Classified image for Riau, Jambi, and South Sumatra using our Combined image stack. Insets focus on the border of Riau and Jambi for (A) Combined, (B) SAR, and (C) Landsat classified images.
Remotesensing 12 01220 g003
Figure 4. Error matrices of the three classifications for Landsat, SAR, and Combined maps. Values in each cell represent the number of testing pixels for each land-cover class.
Figure 4. Error matrices of the three classifications for Landsat, SAR, and Combined maps. Values in each cell represent the number of testing pixels for each land-cover class.
Remotesensing 12 01220 g004
Figure 5. Difference images to compare the agreement of oil palm mapped extents for (A) the Combined-Miettinen2016 map and (B) the Combined-Austin2017 map. Black regions indicate where both maps agree that there is no oil palm, yellow regions indicate where oil palm is present in both maps, blue regions indicate where oil palm is present in the Combined map but not in the reference map (i.e., Miettinen2016 or Austin2017), and pink regions indicate where oil palm is present in the reference map but not in the Combined map.
Figure 5. Difference images to compare the agreement of oil palm mapped extents for (A) the Combined-Miettinen2016 map and (B) the Combined-Austin2017 map. Black regions indicate where both maps agree that there is no oil palm, yellow regions indicate where oil palm is present in both maps, blue regions indicate where oil palm is present in the Combined map but not in the reference map (i.e., Miettinen2016 or Austin2017), and pink regions indicate where oil palm is present in the reference map but not in the Combined map.
Remotesensing 12 01220 g005
Figure 6. Weighted percentage of oil palm land-cover points, classified by their respective maps, obtained from our two disagreement extents (Combined-Miettinen2016 and Combined-Austin2017), which fall under five land-cover classes in Google Earth Pro. Bare ground-crop (BGC); forest shrub mosaic (FSM); forest (FOR); oil palm (OPM); built-up area (URB).
Figure 6. Weighted percentage of oil palm land-cover points, classified by their respective maps, obtained from our two disagreement extents (Combined-Miettinen2016 and Combined-Austin2017), which fall under five land-cover classes in Google Earth Pro. Bare ground-crop (BGC); forest shrub mosaic (FSM); forest (FOR); oil palm (OPM); built-up area (URB).
Remotesensing 12 01220 g006
Table 1. Definition of vegetation indices calculated from Landsat imagery. Near-infrared band (NIR); red band (RED), blue band (BLUE); shortwave infrared band (SWIR).
Table 1. Definition of vegetation indices calculated from Landsat imagery. Near-infrared band (NIR); red band (RED), blue band (BLUE); shortwave infrared band (SWIR).
Vegetation IndicesFormulasCharacteristics
Normalized Difference Vegetation Index (NDVI)NDVI = N I R     R E D N I R + R E D NDVI is widely used as a vegetation metric because it is sensitive to change in greenness and therefore green vegetation biomass.
Enhanced Vegetation Index (EVI) EVI = 2.5     N I R     R E D N I R + 6.0     R E D     7.5     B L U E + 1.0 EVI improves the vegetation signal especially in regions with dense green vegetation. Further, it is little influenced by the soil and atmosphere.
Soil-Adjusted Total Vegetation Index (SATVI)SATVI = S W I R 1     R E D S W I R 1 + R E D + 0.1 ∗(1.1 − S W I R 2 2.0 )SATVI can detect both senescent and green vegetation.
Normalized Difference Tillage Index (NDTI)NDTI = S W I R 1     S W I R 2 S W I R 1 + S W I R 2 NDTI is sensitive to tillage intensity and crop residue cover.
Normalized Difference Water Index (NDWI)NDWI = N I R     S W I R 1 N I R + S W I R 1 NDWI is sensitive to changes in vegetation liquid water content. It is higher for green than dry vegetation.
Table 2. Formulae of indices calculated from SAR imagery.
Table 2. Formulae of indices calculated from SAR imagery.
Radar IndicesFormula
Average (AVE) H H + H V 2
Difference (DIF)HH − HV
Ratio 1 (RAT1) H H H V
Ratio 2 (RAT2) H V H H
Normalized Difference Index (NDI) H H H V H H + H V
NL Index H H H V H H + H V
Angular Second Moment (ASM) i j { p ( i , j ) } 2
Average (AVG) i = 2 2 N g ( i p ( x + y ) ( i ) )
Contrast (CON) n = 0 N g 1 n 2 { i = 1 N g j = 1 N g p ( i , j ) } , | i j | = n
Correlation (COR) i j ( i , j ) p ( i , j ) μ x μ y σ x σ y
Dissimilarity (DIS) n = 1 N g 1 n { i = 1 N g j = 1 N g p ( i , j ) 2 } , | i j | = n
Entropy (ENT) i j p ( i , j ) l o g ( p ( i , j ) )
Inverse Difference Moment (IDM) i j ( 1 1 + ( i j ) 2 ) p ( i , j )
Variance (VAR) i j ( i μ ) 2 p ( i , j )
Table 3. Producer’s (%) and user’s (%) accuracy of Combined, Landsat, and SAR classified maps for all land-cover classes. Producer’s accuracy (PA); user’s accuracy (UA). Values in bold represent the highest producer’s and user’s accuracy for oil palm land-cover and highest overall accuracy.
Table 3. Producer’s (%) and user’s (%) accuracy of Combined, Landsat, and SAR classified maps for all land-cover classes. Producer’s accuracy (PA); user’s accuracy (UA). Values in bold represent the highest producer’s and user’s accuracy for oil palm land-cover and highest overall accuracy.
Land-Cover Class CombinedLandsatSAR
Oil PalmPA84%70%81%
UA90%62%66%
ForestPA76%72%82%
UA94%91%83%
Forest Shrub MosaicPA64%52%27%
UA41%37%24%
Bare Ground CropPA91%86%73%
UA83%80%79%
Built-Up AreaPA97%93%59%
UA97%93%84%
WaterPA98%98%94%
UA99%97%95%
Overall Accuracy84%78%74%

Share and Cite

MDPI and ACS Style

Sarzynski, T.; Giam, X.; Carrasco, L.; Lee, J.S.H. Combining Radar and Optical Imagery to Map Oil Palm Plantations in Sumatra, Indonesia, Using the Google Earth Engine. Remote Sens. 2020, 12, 1220. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12071220

AMA Style

Sarzynski T, Giam X, Carrasco L, Lee JSH. Combining Radar and Optical Imagery to Map Oil Palm Plantations in Sumatra, Indonesia, Using the Google Earth Engine. Remote Sensing. 2020; 12(7):1220. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12071220

Chicago/Turabian Style

Sarzynski, Thuan, Xingli Giam, Luis Carrasco, and Janice Ser Huay Lee. 2020. "Combining Radar and Optical Imagery to Map Oil Palm Plantations in Sumatra, Indonesia, Using the Google Earth Engine" Remote Sensing 12, no. 7: 1220. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12071220

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop