Monitoring tropical forest change using tree canopy cover time series obtained from Sentinel-1 and Sentinel-2 data

ABSTRACT

The most practical method for monitoring forest change over large areas is using remotely sensed data. However, given that current techniques are somewhat weak for monitoring small-scale forest disturbances, achieving accurate monitoring remains challenging, especially in tropical areas where selective and illegal logging occurs frequently. To further improve the ability to monitor forest changes, we estimated tree canopy cover (TCC) using Sentinel-1 and Sentinel-2 data. We developed an approach to monitor forest change on the obtained TCC time series. This approach was applied to monitor forest change in the Bago Mountains of Myanmar from 2017 to 2021. We then completed accuracy assessments and area estimation using reference data obtained from stratified random sampling and unbiased estimators. The final results indicated that: (1) in TCC estimation, Sentinel-1 played a limited role; the red-edge bands of Sentinel-2 achieved slightly different results to the other bands, and superior results were obtained by using all bands; (2) our method successfully mapped forest change with the overall accuracy of 93%. Furthermore, compared with the most widely used and the most recent approaches, our method was better at capturing forest disturbances.

KEYWORDS:

• Deforestation
• forest degradation
• time series analysis
• Sentinel-1
• Sentinel-2
• google earth engine

1. Introduction

Tropical forests account for approximately half of the total forest area and have a significantly higher carbon uptake capacity than other forest types (Pan et al. 2011). However, most tropical forests are located in developing countries (Romijn et al. 2015), where burgeoning economic activities have led to severe deforestation or forest degradation (Pearson et al. 2017). The issue of forest degradation in the tropics is now a focus because it may make up a large portion of carbon emissions from tropical forest loss (Baccini et al. 2017). Tropical countries account for 72% of the total forest area of countries reporting forest degradation (FAO 2020a). Studies on forest degradation have focused on tropical areas in recent years (Bullock, Woodcock, and Olofsson 2020a; Bullock et al. 2020b; Matricardi et al. 2020). However, the focus on tropical forest degradation is uneven at the global level (Dupuis et al. 2020). Most studies examine forest degradation in the Americas but forest degradation remains relatively poorly studied in Southeast Asia, where selective and illegal logging frequently occurs (Miettinen, Stibig, and Achard 2014).

Although field inventories remain the primary method for monitoring forest degradation at the national level (FAO 2020a), collecting these data may be difficult for less economically developed tropical countries, especially those in Southeast Asia with large forest areas and rugged terrain. Satellite remote sensing has unrivaled advantages on macroscopic and dynamic monitoring, thus using satellite remote sensing to monitor forest degradation over large areas is undoubtedly the most practical alternative. However, the subtle and small-scale nature of forest degradation poses a huge challenge for this field. Considerable research has indicated that commonly used Landsat imagery and algorithms have difficulty capturing selective and illegal logging where disturbances frequently occur at a scale smaller than 30 m (Gao et al. 2020; Herold et al. 2011; Shimizu et al. 2017b). This greatly limits the accuracy of forest loss estimates using remote sensing. To overcome this challenge, researchers have made various attempts, for example, some researchers have used spectral mixture analysis (SMA) to obtain the proportion of each endmember in a pixel and then combined the obtained endmember information with commonly used algorithms to successfully capture a certain amount of forest degradation using Landsat imagery (Bullock, Woodcock, and Olofsson 2020a; Chen et al. 2021). However, it is worth noting that most of the existing research on tropical forest degradation has focused on evergreen or semi-evergreen forests (Bullock, Woodcock, and Olofsson 2020a; Bullock et al. 2020b; Langner et al. 2018; Slagter et al. 2023; Vancutsem et al. 2021), and some of these methods may not be as effective in deciduous forests, where monitoring subtle forest degradation could become more complicated due to phenology variations.

The Sentinel-2 (S2) satellite, launched in 2015, provides images with higher temporal and spatial resolution than Landsat, and its red-edge bands (REB) were highly anticipated by the vegetation remote sensing community because of their exceptional sensitivity to chlorophyll content. To date, S2 has been used extensively in land monitoring and most studies have highlighted the positive impact provided by its advantages (Bartold and Kluczek 2023; Korhonen, Hadi, and Rautiainen 2017; Phiri et al. 2020). However, owing to the current amount of data, especially atmospherically corrected data, which is difficult to satisfy the needs of commonly used change detection algorithms, such as LandTrendr (Kennedy, Yang, and Cohen 2010) and CCDC (Zhu and Woodcock 2014), the advantages of S2 have not yet been fully demonstrated in forest change monitoring. Consequently, we need to explore a method for monitoring forest change that can use the available S2 data fully.

Tree canopy cover (TCC) is the estimated percentage of tree cover in a specified pixel and compared with categorical maps, it is more suitable for monitoring forest degradation (Estoque et al. 2021). Some studies have used TCC as an indicator to assess forest gain and loss (Wimberly et al. 2022), and methods for estimating TCC using S2 are becoming abundant (Cilek et al. 2022; Verhegghen et al. 2022). Thus, monitoring forest change using TCC obtained from S2 seems feasible. However, creating TCC time series for monitoring tropical forest change is difficult based on the reported methods. For example, Cilek et al. (2022) obtained reference data using high-resolution commercial imagery, which is not widely accessible for tropical countries with heavy forest loss but underdeveloped economies; whereas TCC obtained by Verhegghen et al. (2022) are categorical data and each category covers a wide range of TCC values, which may lead to an underestimation of forest degradation. In addition, although considerable research indicated that REB was useful in forest monitoring (Astola et al. 2019; Kluczek, Zagajewski, and Zwijacz-Kozica 2023), existing studies do not provide a comparative analysis of the role of REB in TCC assessment. Thus, methods for TCC estimation using S2 data need to be further explored.

The backscatter signal contains information on forest structures and is not susceptible to cloud coverage, and therefore synthetic aperture radar (SAR) has been widely used for forest mapping, above-ground biomass estimation, and disturbance detection (Amitrano et al. 2021). Among SAR data, Sentinel-1 (S1) has great potential because it is the first constellation to make global continuous SAR data freely and publicly available (Torres et al. 2012). However, its usefulness in forest monitoring is considered to be limited owing to its short wavelength (C-band), which limits penetration into the tree canopy (Hirschmugl et al. 2020). Nevertheless, since the launch of Sentinel-1A in April 2014, many studies have used it for forest monitoring and confirmed that it plays a positive role (Hethcoat et al. 2021; Shimizu, Ota, and Mizoue 2019). Different from above-ground biomass, TCC is mainly related to the tree canopy information that can be provided by the C-band (Sinha et al. 2015). Therefore, it is reasonable to attempt to apply S1 to TCC estimation and to assess its usefulness.

The main objective of our study was to develop a method for forest change monitoring that could be applied in underdeveloped tropical areas using S1 and S2 data. The specific objectives were as follows: (1) developing a method for TCC estimation that could meet the requirements of forest change monitoring in the tropics and assessing the potential of S1 and S2 bands and their derived variables for TCC estimation; and (2) developing a workflow for monitoring forest change using TCC time series.

2. Study area

The study area was the Bago Mountains in central Myanmar (Figure 1). The area covers approximately 2.3 million ha, with elevations ranging from 10 to 600 m. The area is dominated by mixed deciduous forests in which the main species are hardwoods with a high economic value (e.g. Tectona grandis and Xylia xylocapa).

Figure 1. Study area. True color composite S2 (20 m) acquired in 2019 was used as a background. National boundaries were based on the World Borders Dataset downloaded from Thematic Mapping (http://thematicmapping.org/).

Selective logging has been conducted since 1856, guided by the Myanmar selection system. However, because of the logging ban (Shimizu et al. 2017a) and ongoing corruption in forest sector (Lim et al. 2017), we were unsure of the state of selective logging in the region after 2016. In addition to legal selective logging, long-standing illegal logging (Saung et al. 2021), extensive shifting cultivation (Ei, Kosaka, and Takeda 2017), and infrastructure development that occurred with the development of Nay Pyi Taw (Lim et al. 2017) have put enormous pressure on the forests.

Notably, a buffer was generated around the study area so that spatial features (see Section 3.3.2) could be calculated for all pixels within the study area.

3. Materials and methods

3.1. Overview

Google Earth Engine (GEE) was used as the data processing platform. Our method had three main parts (Figure 2). The first part was preprocessing the S1 and S2 data to generate S1 and S2 time series according to our needs (Section 3.2). The second part was the generation of the TCC time series using the S1 and S2 time series. In this section, we constructed six different models to analyze the role played by different bands and their derived features in TCC estimation (Section 3.3). The third part focused on mapping forest change using the highest accuracy TCC time series obtained in Section 3.3. We completed the classification through a developed workflow and later followed the good practice steps of Olofsson et al. (2014) to provide unbiased area estimates of forest, non-forest, deforestation, forest degradation, and forest gain in the Bago Mountains (Section 3.4).

Figure 2. Study framework. The S1 time series is composed of RGB visualization images (red: VV, green: VH, blue: VV/VH). The S2 time series is composed of true color composite S2 images.

3.2. Remotely sensed data preprocessing

3.2.1. Sentinel-2 images preprocessing

The S2 Level-2A dataset only provided images from December 2018 for the study area, and therefore we instead acquired all available S2 MSI, Level-1C images from 2016 to 2022. To avoid cloud pollution and to minimize variations in vegetation phenology, only images with cloud cover of less than 30% acquired between 1 November and 31 December were used.

We first implemented atmospheric correction using Satellite Invariant Atmospheric Correction (SIAC) developed by Yin et al. (2019) to reduce the impact of atmospheric effects. The SIAC method can achieve atmospheric corrections for medium-resolution satellites (Landsat-8/Sentinel-2) and it is currently being used in water body related studies (Maciel et al. 2021). Topographic correction is regarded as an indispensable step for mountainous areas (Vanonckelen, Lhermitte, and Van Rompaey 2013). We used the SCS + C model (Soenen, Peddle, and Coburn 2005) for topographic correction. The operation of this model in GEE was developed by Poortinga et al. (2019), and the Digital Elevation Model used was the SRTM V3 product. Next, we masked the cloud and cirrus pixels using the cloud mask band (QA60) and the s2cloudless algorithm (cloud probability: 60), and then resampled the B2, B3, B4, and B8 bands to 20 m using the bicubic interpolation (Runge and Grosse 2020; Zupanc 2017). Finally, we created annual S2 composite images by computing the median values for each year. Any gaps in the S2 time series were filled from averages of the preceding and subsequent years for the same location.

3.2.2. Sentinel-1 images preprocessing

The S1 Ground Range Detected data in the GEE were preprocessed with the S1 Toolbox to perform thermal noise removal, radiometric calibration, and terrain correction. In addition, the refined Lee filter (Lee, Grunes, and de Grandi 1999) was applied to S1 data to reduce the effect of inherent speckle noise. We only selected data from vertical transmit and vertical receive (VV), and vertical transmit and horizontal receive (VH) polarization in descending scenes (Maskell et al. 2021) with the same composite period of S2. Then, we followed the same procedure as S2 to complete the resampling, annual image compositing, and gap filling.

3.3. TCC estimation with random forests

3.3.1. Training and validation data

To obtain reference data, we used Collect Earth, which is a freely available tool developed by the Food and Agriculture Organization (FAO) of the United Nations for monitoring land cover and land use. Using Google Earth technologies, the tool offers access to all available satellite archives of very high resolution (VHR) images (Bey et al. 2016). Because the availability of historical VHR imagery on Google Earth Pro varied over time and space, we constructed a total of 2000 sample plots of 0.04 ha within the study area by simple random sampling only. For each sample plot, we labeled the TCC for all years in which VHR images were available from 1 October to 31 January of the subsequent year. The TCC labeling method was similar to that of Potapov et al. (2019). We created 3 × 3 regular grid points in each sample plot and calculated the proportion of the nine grid points that intersected the tree canopy (Figure S1). Finally, a total of 5010 sample units (1673 in 2016, 717 in 2017, 605 in 2018, 966 in 2019, 315 in 2020, 673 in 2021, and 61 in 2022) were used for training and validation, and we divided the data for each year into training (3508) and validation (1502) data in the ratio of 7:3.

3.3.2. Variables for random forests

We used the random forests (RF) algorithm (Breiman 2001) for regression to obtain TCC maps for each year of the study area. In training the RF model, we used the smileRandomForest function in GEE, where numberOfTrees was set to 500, and the rest of the parameters were kept as the default. To analyze the potential of the S1 and S2 bands and their derived variables for TCC estimation, we constructed six different models in total using 50 variables (Table S1). Three models were basic models, including: (1) a model using only S1-derived indices (here after, s1), (2) a model using B2, B3, B4, B8, B11, B12, and indices derived from these bands in S2 (here after, s2), and (3) a model using REB (i.e. B5, B6, B7, B8A) and indices derived from these bands in S2 (here after, reb). In addition, we constructed three combined models: (4) s1 + s2 model, a model containing all indices of the s1 and s2 models (5) s2 + reb model, a model using all S2-derived indices, and (6) s1 + s2 + reb model, a model using all indices derived from S1 and S2. We used spatial moving windows ranging from 3 × 3 to 15 × 15 pixels to obtain spatial features. To avoid the potential negative impact of redundancy, only the three most important variables in the s1, s2 and reb models (before adding spatial features) were selected for spatial statistics.

3.3.3. Accuracy assessment of TCC estimation

To assess the accuracies of the models in TCC estimation, we calculated the root mean squared error (RMSE) and coefficient of determination (R²) for each model using validation data (Bennett et al. 2013). The RMSE measures the closeness of the prediction to the observation, with lower values indicating higher accuracy, while the R² reflects how much the independent variable explains the variation, with higher values indicating better results.

3.4. Classification of land cover/forest change types

3.4.1. Land cover and forest change classes

Given that there are various definitions of forest degradation but the definitions of deforestation are generally consistent, in this study, forest degradation was defined as a decline in TCC owing to human or natural causes other than deforestation. Furthermore, we classified forest degradation into two types depending on whether recovery signs were present. The recovery threshold was determined as 70% of mean TCC value before disturbance (Shimizu et al. 2022). The seven classes that we constructed are shown below:

• Forest: A pixel with TCC greater than or equal to 40% and without disturbance throughout our study period, where 40% is conservatively defined by referring to the FAO’s definition of closed forest in Myanmar (FAO 2020b).

• Non-forest/Other: Not meeting the definition of other classes, including non-forest area (e.g. farmland, water body, infrastructure), change in other land covers, and forest disturbance after forest gain or recovery.

• Deforestation (Def): A pixel experiencing permanent or long-term transformation from forest to non-forest area, excluding pixels that still had recovery potential after clear-cutting.

• Degradation with recovery (Deg_rec): A pixel experiencing forest degradation and emerging with identifiable nascent tree canopy in our study period.

• Degradation without recovery (Deg_no_rec): A pixel experiencing forest degradation but not meeting the definition of Deg_rec.

• Forest gain: A pixel that was not a forest at first but developed into forest and remained forest at the end of our study period.

• Buffer: 1-pixel buffer around the Def class. We defined this class because a study by Olofsson et al. (2020) showed that a spatial buffer stratum can mitigate the uncertainty of omission errors in small change classes that occur within a large class (e.g. forests).

3.4.2. Classification workflow

We constructed a threshold-based classification workflow to classify land cover/forest change types for each year from 2017 to 2021 (Figure 3). The variables used in the workflow were derived from the TCC time series, while the thresholds we used, except the threshold used to distinguish disturbance (Losthr) and deforestation (Defthr), were decided based on the class definitions. To determine Losthr and Defthr, we investigated 1344 pixels (Forest: 558, Def: 220, Deg_rec: 178, Deg_no_rec: 388), where class could be identified by VHR imagery in Google Earth Pro and the years of disturbance varied. We determined Losthr and Defthr to be 0.3 (Figure S2) and 7.5% (Figure S3), respectively, with reference to the overall accuracy, commission error, and omission error when the two variables were changed. After obtaining the classification maps for each year, we composited the annual classification maps based on the principles that forest disturbance over other types, deforestation over forest degradation, and earlier disturbances over later disturbances, as well as the essential logical relationships between classes, to obtain the final map of classification results.

Figure 3. Classification workflow. (Focal: TCC value of the focus year; Pre: TCC value of the year before the focus year; Pos: TCC value of the year after the focus year; Max: maximum value of TCC time series; MinBefore: minimum value of TCC before the focus year; MaxBefore: maximum value of TCC before the focus year; MinAfter: minimum value of TCC after the focus year; MaxAfter: maximum value of TCC after the focus year; MeanAfter: mean value of TCC after the focus year; LossRate: TCC loss rate = (Pre-Focal)/Pre; Forthr: forest threshold = 40%; Recthr: recovery threshold = 70%MeanBefore; Prethr: threshold for suppressing incorrect disturbance = 95%Pre).

3.4.3. Accuracy assessment and area estimation

We completed accuracy assessment and area estimation for all classification maps by stratified random sampling as recommended by Olofsson et al. (2014). For the final classification map, the sampled population was all of the 20 m × 20 m pixels covering the study area, and the sample size was determined by the formula provided by Cochran (1977, Equation (5.25)). We set the target standard error of overall accuracy (OA) to 0.01 and determined the final sample size to be 1655. To allocate sample units by strata, we ensured a minimum sample size of 30 (Olofsson et al. 2020) while allocating them in close proportion to the area of each stratum (Table 1).

Table 1. Strata used for accuracy assessment.

Stratum	Forest	Non-forest	Deg_rec	Deg_no_rec	Def	Forest gain	Buffer	Total
Map area (km^2)	13974.478	6247.099	1346.689	804.295	61.539	463.452	204.004	23101.556
Strata weights	0.605	0.270	0.058	0.035	0.003	0.020	0.009	1.000
Sample allocation	975	435	95	55	30	35	30	1655

To further verify the effectiveness of the method used in this study, we compared the TCC-based results with the Global Forest Change (GFC; Hansen et al. 2013) dataset and with the results using CCDC-SMA (Chen et al. 2021). The GFC is one of the most widely used global forest change dataset, and the CCDC-SMA is the most recent approach to monitoring forest degradation developed by combining CCDC and SMA. Because the GFC and CCDC-SMA-based results did not directly provide information about disturbance types and forest gain, we only analyzed the ability of the three methods to detect forest disturbance. We obtained the data for GFC and CCDC-SMA-based results on GEE. Forest disturbances detected by GFC outside the period 2017–2021 were reclassified as no disturbance. The CCDC-SMA-based results were obtained using the Tropic Collection 2 (https://github.com/shijuanchen/forest_degradation_georgia), where the required forest mask was constructed from the 2016 TCC results (forest threshold of 40% and nearest neighbor resampling to 30 m), and the required parameters were kept as the default. TCC-based results were reclassified into two classes: disturbance (including Deg_rec, Deg_no_rec, and Def) and no disturbance (all remaining classes). On the basis of the reclassification results, we obtained 1500 pixel locations (no disturbance: 1282, disturbance: 218) at 30 m scale using the stratifiedSample function.

VHR images in Google Earth Pro, PlanetScope images, and S2 (10 m) images were used to label the class and associated year for each sample point for the accuracy assessment. After labeling, we estimated the accuracy, area, and confidence interval of each stratum for the final classification map, and the accuracy of the three disturbance/no disturbance maps. Because the disturbance year may vary between maps owing to the effects of the detection method, the composite period of images, and the availability of high-resolution images, we allowed for a 1-year temporal difference for the correct disturbance detection in assessing accuracy of each disturbance/no disturbance map (Shimizu, Ota, and Mizoue 2020).

4. Results

4.1. TCC estimation

Of the six models, the s1 model had the worst results (R² and RMSE were 0.36 and 35.06, respectively), with the remaining two basic models showing similar accuracy. The combined models all showed some improvement in accuracy compared with the basic models, but the improvement was limited. The s1 + s2 + reb model performed best (R² and RMSE were 0.9 and 14.24, respectively) (Table 2).

Table 2. Accuracy assessment results of different models.

Accuracy Statistic	s1	s2	reb	s1 + s2	s2 + reb	s1 + s2 + reb
R^2	0.36	0.88	0.87	0.89	0.89	0.90
RMSE	35.06	15.44	15.60	14.86	14.67	14.24

4.2. Accuracy assessment and area estimation using the final classification map

The final classification and year of change maps were obtained using the classification workflow (Figure 4). The OA of our maps was approximately 93%, with high accuracy for both forest and non-forest classes and relatively low accuracy for change classes (Table 3). Although the accuracy of the forests was high, the user’s accuracy (UA) of the forests was significantly higher than the producer’s accuracy (PA). The margin of error (defined as the half width of the 95% confidence interval divided by the area estimate) for area estimates was less than 100% for all classes but was large for the change classes, where Deg_rec, Deg_no_rec, Def, and Forest gain were 15%, 18%, 59%, and 24%, respectively. For comparison with other studies, we also merged degradation with and without recovery into the degradation class, which ultimately yielded a PA and UA of 86% and 71%, respectively, for the degradation class (Table S2).

Figure 4. Map of land cover/forest change types from 2017 to 2021 for the entire study area with detailed information (Sentinel-2 in 2016 and 2021 are true color composite images).

Table 3. Confusion matrix consisting of proportions, unbiased area estimates, and accuracies for all classes (area std. error: standard error of area estimates (km²)).

Reference
Map strata	Forest	Non-forest	Deg_rec	Deg_no_rec	Def	Forest gain	Total	UA	PA
Forest	0.594	0.002	0.002	0.002	0.000	0.003	0.605	0.983	0.945
Non-forest	0.009	0.256	0.001	0.002	0.000	0.003	0.270	0.945	0.947
Deg_rec	0.013	0.004	0.039	0.001	0.000	0.002	0.058	0.663	0.873
Deg_no_rec	0.006	0.002	0.003	0.024	0.000	0.000	0.035	0.691	0.729
Def	0.000	0.001	0.000	0.001	0.001	0.000	0.003	0.367	0.624
Forest gain	0.003	0.002	0.000	0.000	0.000	0.014	0.020	0.714	0.632
Buffer	0.002	0.003	0.000	0.003	0.001	0.000	0.009
Total	0.629	0.270	0.044	0.033	0.002	0.023	1.000	OA
Area Estimates (km^2)	14526.064	6230.168	1023.254	762.072	36.165	523.834
1.96 × (area std. error)	221.012	178.063	153.620	134.626	21.437	123.498		0.928
Margin of error	0.015	0.029	0.150	0.177	0.593	0.236

The accuracy comparison results for the three maps are shown in Figure 5. Benefiting from excellent disturbance detection ability (PA: 87%) and low false detection (UA: 67%), the TCC-based map achieved the highest OA of 96%. The GFC had the second highest OA (89%), mainly because of its superior classification of large strata (no disturbance), but it had the poorest detection ability for disturbances, with a PA of only 14%. The PA (61%) of CCDC-SMA was significantly higher than that of GFC, however, its OA (86%) was slightly lower than that of GFC owing to too many false detections in disturbance (UA: 44%).

Figure 5. Accuracy comparison of the GFC, CCDC-SMA, and TCC-based maps: (a) OA; and (b) disturbance class. The error bars represent the 95% confidence interval.

The area of forest disturbance caused by deforestation was tiny, with a decrease followed by an increase from 2017 to 2021, and the deforestation area in 2017 was significantly higher than that in other years. Among the forest degradation types, degradation without recovery showed an increasing trend year by year (in addition to 2018). Accordingly, degradation with recovery was dominant in forest degradation, except in 2021. The area of forest gain was highest and significantly greater in 2017 than in other years, and this phenomenon was probably related to the recovery from forest disturbance that occurred in 2016. Since 2017, forest gain has decreased rapidly and shows a trend towards stabilization (Figure 6).

Figure 6. Annual area of forest change using pixel counts from 2017 to 2021.

5. Discussion

5.1. TCC estimation using Sentinel-1 and Sentinel-2

In TCC estimation, the accuracy of using S1 variables alone was low. This may be related to the properties of S1 and the method for estimating TCC. Nicolau et al. (2021) indicated that obtaining high-accuracy classification results in the tropics using the tree classifier is challenging because some land cover classes (e.g. farmland, secondary vegetation) have similar backscatter responses. Meanwhile, because C-band is susceptible to rapid saturation in dense tropical forests (Meyer 2019), the acquisition of continuous TCC data may be more complex than classification with specific classes. Although using S1 variables alone could not accurately estimate TCC, the addition of S1 variables helped to increase the accuracy of TCC estimation and final classification. In our study, some pixels had premature leaf discoloration. For these pixels, the s1 + s2 + reb model tended to give higher TCC values than the s2 + reb model (Figure S4), which ultimately had a positive effect on reducing commission errors of disturbance classes (especially for Def with a small area). This result may be related to the lack of sensitivity of the SAR data to leaf color changes. However, because the difference between the s1 + s2 + reb model and the s2 + reb model was not obvious (Table 2; Figure S5), future studies can consider using only S2 data for TCC estimation to reduce data processing time and cost.

Both the s2 model and the reb model could estimate TCC values well, but there were clear differences in their results for typical land covers. For example, the TCC value obtained by the s2 model reached approximately 80% one year after a clear-cutting operation, which is highly unlikely to happen in reality. Conversely, the reb model gave more realistic results (Figure S6), which may be related to the fact that REB have a higher association with plant physiological traits (Korhonen, Hadi, and Rautiainen 2017; Zarco-Tejada et al. 2018). However, probably because of the original spatial resolution, the reb model was not ideal when rivers and trails existed in closed forests compared with the s2 model, which caused the reb model to overestimate the TCC values of the corresponding pixels (Figure S5(b)). Both the s2 model and the reb model had certain shortcomings. In contrast, the combined model gave better results (Figure S5). Therefore, we suggest that all available bands of s2 should be fully used in estimating TCC if conditions permit.

5.2. Forest change detection using TCC time series and threshold classification

We monitored forest changes using TCC time series and threshold classification. The OA of our classification map was 93%, with a PA and UA of 86% and 71% for the forest degradation class, and 62% and 37% for the deforestation class, respectively. The study conducted by Vancutsem et al. (2021) is quite similar to ours, and the study area also covers the Bago Mountains, which could have allowed for direct comparisons with our study. Regrettably, this study classified deciduous forests as other land cover and could not provide information on forest change for the corresponding class. Kyaw et al. (2020) investigated forest loss and forest degradation in parts of our study area from 2000 to 2017, and indicated that forest degradation paralleled forest loss in the relevant regions, however, this result differs significantly from our study (Figure 6). We suggest that there are two main reasons for this difference. Firstly, Kyaw et al. (2020) defined the main disturbance types (e.g. timber plantation, shifting cultivation) as forest loss, while in our study the relevant areas were defined as forest degradation; secondly, there was extensive dam construction in the early twenty-first century, but during the study period of our study, deforestation was mainly dominated by small-scale infrastructure construction. Therefore, in contrast to the view of Kyaw et al. (2020) our results suggest that the primary forest disturbance in the Bago Mountains is forest degradation rather than deforestation. Compared with similar studies in other regions (Bullock, Woodcock, and Olofsson 2020a; Bullock et al. 2020b; Chen et al. 2021), our PA of forest degradation was higher, which indicated to some extent that our method was more capable of capturing forest degradation. However, the accuracy of our deforestation class was slightly lower than other studies, but it should be noted that there was a difference in the definition of deforestation and forest degradation between our study and previous studies. As mentioned earlier, in our study, deforestation excluded pixels that still had recovery potential after clear-cutting, which was more consistent with the FAO standard (FAO 2012). Meanwhile, previous studies used indicators (e.g. intercept, regression coefficients) obtained from relevant models to detect deforestation and forest degradation (Bullock, Woodcock, and Olofsson 2020a; Chen et al. 2021). However, for practicality, these indicators should be the same as those used by countries and regions to define forest degradation and deforestation. In this study, unlike previous studies, TCC was first estimated and then provided to the workflow; TCC is an indicator used in FAO’s definition of forests and is considered practically useful.

The comparison results for the three maps showed an obvious difference between the three methods for disturbance detection capability, which is demonstrated visually in Figure 7. The TCC-based approach was great for disturbance detection at various scales (ranging from approximately 0.05 ha to over 50 ha). The CCDC-SMA was the second best, and it could detect large-scale disturbance well, while also capturing some sub-pixel and relatively small-scale disturbances. The GFC was similar to CCDC-SMA in detecting large-scale disturbance but essentially could not detect relatively small-scale disturbance. However, with the improvement in disturbance detection capability, false detections also markedly increased. The false detections in the TCC-based model were mainly attributed to premature leaf discoloration and abscission (Figure 7(d)), causing considerable changes in reflectance in different years. The CCDC-SMA model was also affected by phenological anomalies, but there was a large number of other false detections with unclear causes, which could be related to the greater challenges of model fitting in deciduous forests (Chen et al. 2021) and the lack of parameter optimization. The GFC model achieved the lowest number of false detections and the highest reliability of detected disturbances.

Figure 7. Differences comparison between TCC-based, CCDC-SMA, and GFC: (a) small-scale forest disturbances (17.817°N, 96.458°E); (b) frequently disturbed area (18.448°N, 95.849°E); (c) large-scale logging (18.238°N, 95.903°E); (d) region of phenological anomalies (20.012°N, 95.857°E). True color composite VHR images from Google Earth Pro are used as a background.

Among the three methods, the acquisition of GFC was the most convenient, but it was weakest at detecting small-scale disturbance and severely underestimated the actual forest loss. The CCDC-SMA could detect sub-pixel forest disturbances with only a forest mask provided by the user, but the results obtained by this method notably had a high number of false detections with the default parameters. Therefore, parameter optimization and terrain correction (Miettinen, Stibig, and Achard 2014) may be necessary if this method is used to obtain reliable results in Southeast Asia. The method used in this study was relatively complicated because it required high quality TCC reference data or TCC time series, but the accuracy of the classification map obtained by this method was high. Our method was completed using GEE, and the parameters and thresholds used in the classification process could be flexibly modified according to actual needs, with the possibility of further extension and application in other study areas.

Although our method achieved high overall accuracy and was good at detecting small-scale disturbances, it still had some limitations. First, we did not fully address the problem of high commission errors in forest degradation classes caused by phenological anomalies. While we tried to suppress the generation of commission errors by introducing additional parameters (Prethr) in the classification process, the effect was limited. Therefore, future research could focus on the image composite stage (i.e. adjusting the image composite period or exploring more suitable image composite methods) or performing post-processing. Second, the accuracy of our deforestation assessment was low. There were two main reasons for this: (a) in the initial years of the time series, the TCC of some farmland pixels was overestimated, which resulted in the corresponding pixels being misclassified as deforestation; and (b) in the terminal years of the time series, the TCC of some Def and Deg_no_rec pixels was similar, which resulted in misclassification between these two classes. Both of these were related to TCC estimation in some way, and thus future research will need to explore methods to improve the accuracy of TCC estimation further.

6. Conclusions

We completed TCC estimation at 20 m spatial resolution using S1 and S2 data, and then we mapped forest change in the Bago Mountains using the obtained TCC time series and the classification workflow. S1 data could not accomplish TCC estimation alone, but its addition did improve the accuracy; using the REB or the other bands of S2 alone could achieve good results, but results using all bands tended to be better. The OA of our forest change maps was 93%, and the accuracy for most land classes was similar to other studies; however, our method was much better at capturing small-scale forest disturbance. This result indicates that monitoring forest change using TCC time series and our classification workflow is feasible and meaningful.

Acknowledgments

We thank Leonie Seabrook, PhD, from Edanz (https://jp.edanz.com/ac) for English editing.

We thank anonymous reviewers for their constructive comments.

Data availability

The data supporting this study will be provided by the corresponding author upon reasonable request.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This study was funded by JSPS KAKENHI [grant number JP19H04339] and Grant for Environmental Research Projects from the Sumitomo Foundation, Japan.