Drivers of Forest Loss in a Megadiverse Hotspot on the Paciﬁc Coast of Colombia

: Tropical forests are disappearing at unprecedented rates, but the drivers behind this transformation are not always clear. This limits the decision-making processes and the e ﬀ ectiveness of forest management policies. In this paper, we address the extent and drivers of deforestation of the Choco biodiversity hotspot, which has not received much scientiﬁc attention despite its high levels of plant diversity and endemism. The climate is characterized by persistent cloud cover which is a challenge for land cover mapping from optical satellite imagery. By using Google Earth Engine to select pixels with minimal cloud content and applying a random forest classiﬁer to Landsat and Sentinel data, we produced a wall-to-wall land cover map, enabling a diagnosis of the status and drivers of forest loss in the region. Analyses of these new maps together with information from illicit crops and alluvial mining uncovered the pressure over intact forests. According to Global Forest Change (GFC) data, 2324 km 2 were deforested in this area from 2001 to 2018, reaching a maximum in 2016 and 2017. We found that 68% of the area is covered by broadleaf forests (67,473 km 2 ) and 15% by shrublands (14,483 km 2 ), the latter with enormous potential to promote restoration projects. This paper provides a new insight into the conservation of this exceptional forest with a discussion of the drivers of forest loss, where illicit crops and alluvial mining were found to be responsible for 60% of forest loss. inﬂuence accuracy estimates or the overall trends. Forest cover is the most abundant class and is the one with the highest likelihood to express strong spatial autocorrelation. We calculated the distance between training and validation points of the most important class, the forest class; from 1969 training points 102 are within 100 m of validation points, 542 are within 500 m of validation points and 820 are within 1000 m of validation points. The rest, 1149 training points have a distance above 1 km from the validation points. In absence of a rigorous assessment, we advise caution about the potential inﬂuence of autocorrelation in some of the results obtained here.


Introduction
Tropical forests are known as significant reservoirs of carbon and play an important role to reduce the negative impacts of climate change as indicated by the project for Reducing Emissions from Deforestation and forest Degradation (REDD+). Tropical forests also sustain a large number of species and therefore are essential for biodiversity conservation [1][2][3]. Despite the importance of the Pacific forests of northern South America in terms of carbon storage [4] and outstanding biodiversity [2,5], important knowledge gaps remain about the status and trends of forest loss in this region, mainly caused by the absence of wall-to-wall forest monitoring programs [6][7][8]. Preserving natural ecosystems requires a timely identification of priority forest areas for conservation and restoration that can provide ecosystems services and that at the same time, are undergoing great pressure [9]. In particular, areas associated with armed conflicts have been related to forest fragmentation due to legal and illegal land uses, including perennial crops, cattle ranching, mining, selective logging and illicit crops [10]. Therefore, information on the status and change of natural ecosystems is essential for implementing conservation policies in light of emerging challenges associated with sociopolitical events.
The future conservation of sensitive ecosystems in Colombia is a source of concern, given the latest sociopolitical events related to the signature of peace agreements with illegal armed groups [8,11]. The government estimates that because of the peace agreement with the Revolutionary Armed Forces of Colombia (FARC), a major Colombian rebel group, agriculture will increase from 25,000 ha in the last three years to more than 1 million ha in the near future [12].
There are also unexpected results associated with the peace process because of changes in territorial control [13] by illegal groups which frequently finance themselves with revenues from illicit crops and alluvial mining, creating negative impacts in the form of deforestation and damage to ecosystems [14]. Unfortunately, official land cover maps are not regularly updated in the country, constraining the information available to make informed decisions.
Forest degradation constitutes an increasing pressure for conservation of biodiversity in Colombia's Pacific Coast and their protected areas [15], due to several drivers of forest loss [4,16] that include land tenure, illegality, and law enforcement [17]. However, the magnitude of the transformation associated with each driver is unknown and therefore quantifying such relations is imperative. This paper aims to: (i) generate a spatially continuous land cover map in an area of high cloud cover persistence; (ii) combine land cover and Global Forest Change (GFC) data [18] to identify the influence of land cover in forest loss, and iii) analyze forest loss in terms of alluvial mining and illicit crops.

Study Area
About half of the world's tropical forests are located in South America. The Amazon Basin forest, Brazilian Atlantic forest and Pacific forest (Ecuador, Colombia, Panama and Costa Rica) are the largest continuous ecoregions of forest in the continent [19]. The Pacific ecoregion is located in the intertropical convergence zone (ITCZ) with altitude gradients varying from the coastal tropical forests to the high Andean páramos located 3000 masl. This area ( Figure 1) is a biodiversity hotspot due to high levels of biological endemism [2]. The area is characterized by persistent cloud cover and is considered as one with the largest annual precipitation globally.
The main rivers are the Patía, San Juan, Baudó, and Mira, flowing towards the Pacific Ocean; and the Atrato River flowing towards the Caribbean Sea. The climate of this tropical rain forest has a short or absent dry season [19] where 55% of the study area has precipitation values of over 5000 mm/year [20]. Under these conditions, the main continental ecosystems are "terra firme" forests as well as coastal vegetation forming mangroves [21,22], coastal marsh ecosystems [23], and wetland forests [4]. There are also transformed areas, including large extensions of industrial plantations to the north and to the south of the study area, as well as other land uses including cattle grazing, growing coca crops, commercial logging and alluvial gold mining.
The 99,852 km 2 of the study area ( Figure 1b) intersects four departments: Chocó, Valle del Cauca, Cauca and Nariño. These departments recognized collective land rights for Afro-Colombian (52,225 km 2 ), and indigenous groups (17,973 km 2 ) through a land reform that took place in the 1990s. According to the law, these collective lands are imprescriptible, inalienable and indefensible [24].

Data and Pre-Processing
Monitoring land cover/land use and determining drivers of biodiversity loss using optical remote sensing is challenging in the context of humid/pluvial forests due to the high or permanent cloud cover [13,14]. To overcome this limitation, we used both optical and synthetic aperture radar (SAR) data for the land cover classification that included information from Landsat-8, Sentinel-2 and Sentinel-1. The GFC database [15] was associated with the updated land cover map in order to analyze the fate of deforested pixels ( Figure 2).
Modern Earth Engine computing capabilities enable us to select the best land cover observations from Landsat and Sentinel satellites, minimizing invalid observations due to clouds' persistence. The time series of optical data helped to identify observations free from clouds during a given period, while the use of active SAR datasets helped to overcome such limitations given their ability to penetrate clouds and haze and their lower sensitivity to atmospheric conditions [25][26][27]. We built optical and radar mosaics for the study area (the size of Portugal) and then performed land cover classification and validation through the implementation of a random forest algorithm [28] using the R programming environment [29].
Data from optical and SAR sensors were fused for the forest and land cover mapping [30]. Given the frequent cloud cover in the area, we included four optical sensors (Landsat 8 OLI, Landsat 7 ETM+ and Sentinel 2-A/B MSI). Several years of optical data acquisitions (2014-2018) were required to obtain a cloud and shadow free composite [31]. Radar wavelengths from Sentinel 1 C-band was also included to complement optical composites. We used the massive computational capabilities of Google Earth Engine (GEE) to access and process multiannual composites [32]. A shuttle radar topography mission (SRTM) layer was also added to the spectral data because it increases classification accuracy [31]. and Sentinel 2-A/B MSI). Several years of optical data acquisitions (2014-2018) were required to obtain a cloud and shadow free composite [31]. Radar wavelengths from Sentinel 1 C-band was also included to complement optical composites. We used the massive computational capabilities of Google Earth Engine (GEE) to access and process multiannual composites [32]. A shuttle radar topography mission (SRTM) layer was also added to the spectral data because it increases classification accuracy [31].  After analyzing the presence of data gaps and speckle, the ascending orbit with VV polarization was selected as an explanatory image band; cross-polarization (HV or VH) is known to result in weaker backscatter [33]. All available images from 2014-2018 were selected in order to obtain the median VV backscatter coefficient (σ • ) in dB. Time series of SAR data are commonly used to minimize speckle and remove extreme backscattering values that usually constitute noise. The C-band is sensitive to the scattering elements and their electromagnetic characteristics, with forests producing high backscatter (scattering radiation towards the SAR sensor) while grassland or water bodies scatter more incident radiation away from the SAR sensor. Sentinel-1 data were processed using the GEE platform based on the COPERNICUS/S1_GRD database. Available observations were filtered to select ascending passes with VV polarization, and the median backscattering per pixel was selected as the radar metric for classification. We assigned to each pixel the median backscattering during the four-year time period as an input for the classification process. The median has been proven to be a more robust metric against extreme values than the mean [34][35][36].

Landsat 7-8 and Sentinel-2
Clouds affect the usefulness of optical pixels and a time series is required to select the best observations by filtering clouds and cloud shadows [37]. The GEE cloudScore principles were applied to Landsat and Sentinel data [38], including the analysis of solar geometry to detect shadows; once these low quality pixels were filtered a median reducer was used for compositing. Top of atmosphere (TOA) reflectances of the Red, NIR, SWIR1 and SWIR2 bands were selected with a 30 m output resolution in order to match the lower resolution of Landsat data (Table 1). Spectral integration of these multiple satellite systems is feasible and discussed in the scientific literature [38][39][40].

Training and Validation Data
Training and validation data were produced through a random stratified design. Strata consisted of the classes considered in the corine land cover (CLC) classification maps for Colombia in 2011 [42]. The number of points sampled within each stratum was proportional to the extent of the stratum in the study area. This resulted in the collection of 9708 interpreted points associated with eight land cover classes. Observed and interpreted points were randomly divided into two groups of 6796 points (70%) for training and 2912 points (30%) for validation ( Table 2). Table 2. Visually interpreted points used to train the random forest classifier and validate the resulting land cover map using the international geosphere biosphere programme (IGBP) classes.

Class Points Train Validate
Bare soil  895  621  274  Cropland  845  610  235  Grassland  429  297  132  Broadleaf forest  2816  1969  847  Shrubland  2552  1799  753  Wetland forest  786  541  245  Wetland grassland  636  429  207  Water Bodies  749  530  219 We assigned a land cover class to each sampled point through field observations using a hand-held GPS GARMIN Map 78 from 2016 to 2018 and through the interpretation of high-resolution images and the official land cover map. Fieldwork was concentrated near the Atrato River due to security restrictions, emphasizing the proper identification and image interpretation of areas affected by alluvial mining, shrublands, and industrial crops. Only large-scale plantations such as banana or oil palm were interpreted as the cropland class. Small-scale agriculture or dispersed crops were not differentiated due to the ambiguous/mixed spectral response. Two additional classes were included from official maps: urban-build-up and páramos. Páramos are high Andean ecosystems located in Peru, Ecuador, Colombia and Venezuela. Colombia has an official delimitation of these ecosystems given their importance in terms of water supply and their high levels of biological endemism [43][44][45].

Classification and Accuracy Assessment
The random forest algorithm was used to classify remote sensing data (Red, NIR, SWIR1, SWIR2, Sentinel-1 and SRTM) into eight land cover classes (Table 2). Explanatory variables were obtained from Remote Sens. 2020, 12, 1235 6 of 16 Landsat 8, Sentinel-1, Sentinel-2 and the SRTM. Random forest classification is a supervised algorithm based on ensemble learning.
Ensemble learning combines the outcomes of different iterations of the same algorithm to improve prediction power. Random forests are comprised by multiple decision trees that are created using a subset of random samples from the training set based on attributes of the dataset of interest. Each decision tree produces a result for each pixel that is considered as a vote. The final prediction for each pixel corresponds to the most voted class among all the decision trees [28,46]. For this application, we used the Random Forests package [47] implemented in the R computing environment. The model was parameterized to run 500 trees using a random selection of variables in each split equal to the square root of the total number of predictors. The accuracy of the resulting land cover map was assessed by using the validation data and the statistics of the confusion matrix.
The land cover classification legend corresponded to the IGBP classes [48]. Finally, this land cover map was overlaid with the forest loss pixels identified by the GFC [18] in order to determine the current land cover condition of deforested pixels from 2001 to 2018.

Ancillary Data
We used estimates of forest loss from GFC. This dataset is available in GEE as Hansen Global Forest Change v1.6. All pixels labeled as forest loss from 2001 to 2018 were selected. To increase our understanding of the role of the drivers of forest loss in the magnitude of the deforestation we included two of the most relevant drivers: illicit crops and alluvial mining, which are highly related to the economy of local communities [49]. The World Drug Report of the United Nations Office for Drugs and Crime (UNODC) [50,51] estimated that 69% of the global coca bush crops are located in Colombia, with an extension of 1461 km 2 in 2016 from which 577 km 2 (39.5%) are located in the Pacific coast.
Spatial explicit layers at 1-km 2 grid available from UNODC reports were combined with deforested pixels using zonal statistics. Appendix A provide web links that point to layers showing areas under the influence of alluvial mining and illicit crops; these drivers of forest loss were treated as zones within the map algebra syntaxes in order to calculate statistics of forest loss. In addition, each class from the land cover map derived from Section 2.3 was treated as a zone in order to associate each forest loss pixel to the current land cover class.

Results
This study used all the available images of Landsat-8 and Sentinel-2 from January 2014 to December 2018 for the study area. The classified land cover map had an overall accuracy of 82% ( Table 3). The accuracy results represent a strong agreement despite the limited number of observable pixels available and the large temporal compositing method required in these humid-pluvial forests.
The main land cover is broadleaf forest, covering 67,473 km 2 (Table 4), equivalent to 68% of the area, followed by shrubland class, covering 14,483 km 2 ,equivalent to 15% of the area (Figure 3). These two classes cover a sum of 81,956 km 2 . Wetland broadleaf forests and wetland grasslands are also significantly represented classes with 7282 km 2 and 4355 km 2 respectively ( Table 4). The largest wetland is associated with the Atrato River flowing into the Caribbean Sea, followed by the large coastal mangroves along the Pacific coast. Grasslands are mainly located to the north, dedicated to cattle grazing, while shrubland includes natural regeneration of abandoned lands, shifting agriculture [52] and illicit crops [53]. Grasslands and shrublands were the most challenging classes to classify due to the different stages of natural regeneration, resulting in low accuracy values (Users' Accuracy 74% and 73% respectively) ( Table 4). Shifting agriculture, natural regeneration and coca crops were included in the shrubland class since these highly fragmented mosaics were not spectrally distinct. Coca crops are highly fragmented and mixed with shifting agriculture and pastures; based on UNODC [54] more than 60% of coca is grown in small plots (average 0.96 ha). Additionally, large-scale producers scatter illicit crops in the landscape to avoid aerial/terrestrial reconnaissance, resulting in a mosaic of natural forests with small and dispersed cultivated fragments. Bare soils included classes such as degraded lands due to poor agricultural practices, exposed gravel/sediments from riverbeds and abandoned mining sites [55]. Detection of mining sites is challenging because their reflectance is similar to other classes such as bare soils (barren) occurring naturally, sediment bar depositions along the rivers, or bare soil from agricultural activities before planting.
Industrial crops are banana and oil palm plantations in the north and south of the study area, respectively. These industrial croplands are related to the construction of channels that drain the excess of surface water [56]. This land use represents less than 1% of the total area with an extent of 854 km 2 . Another important agricultural land use is cattle grazing on grasslands, with an extent of 1177 km 2 established by clearing natural forests, which is a common practice in the neotropics [57]. A major limitation for cattle are the high levels of precipitation. As expected, water bodies have an important extent (2121 km 2 ) in this pluvial environment; rivers and wetlands are fed by watersheds that reach the Páramos class (671 km 2 ). The total deforested pixels from 2001 to 2018 was equivalent to an area of 2324 km 2 . Each deforested pixel was related to the current land use in order to track the fate of deforested pixels over this time.
The largest extent of deforestation (annual and total) from 2001 to 2018 was related to the shrubland class with a total of 1915 km 2 (Table 5), followed by wetlands (206 km 2 ), grasslands (108 km 2 ) and cropland (56 km 2 ). From 261 km 2 pixels as deforested pixels in 2017, 225 km 2 are currently under the shrubland class (Table 5). Finally, deforested areas connected to grasslands were higher than deforested pixels connected to cropland, indicating that deforested pixels are more likely to become grasslands than industrial plantations. Based on UNODC delimitation, areas under the influence of alluvial mining, illicit crops or both, include 25,609 km 2 of the study area. Forest loss where mining occurs sum to 140 km 2 , while forest loss where illicit crops occur sum to 1181 km 2 , and forest loss where both drivers occur sum to 77 km 2 . Both drivers were present during the 2001-2018 period. Similar patterns of the drivers of forest loss ( Figure 4) are partially due to the fact that areas of illicit crops and alluvial mining overlap, indicating their coexistence in this biodiversity hotspot. A major change in this pattern was found for 2018 when forest loss increase was associated mainly with illicit crops.
Remote Sens. 2020, 12, x FOR PEER REVIEW 9 of 16 clearly interpreted from optical composites or radar images (Figure 3), where more than 100 km of riparian vegetation was removed from both sides of the Quito River. The shrubland class occupied 14,483 km 2 , (see Table 4), which includes abandoned crops, abandoned pastures, and slash and burn plots. A major use of this class also includes illicit coca crops, which reached a maximum of 656 km 2 in 2017 [58].   From the 282 km 2 deforested in 2016 in the Pacific region ( Figure 5), 70% of mining activity was related to alluvial mining along the Atrato River, where heavy-duty machinery (excavators and dredges) is used for vegetation removal and soil extraction. Degradation due to alluvial mining is clearly interpreted from optical composites or radar images (Figure 3), where more than 100 km of riparian vegetation was removed from both sides of the Quito River. The shrubland class occupied 14,483 km 2 , (see Table 4), which includes abandoned crops, abandoned pastures, and slash and burn plots. A major use of this class also includes illicit coca crops, which reached a maximum of 656 km 2 in 2017 [58].

Discussion
Official land cover maps are not regularly updated in Colombia, making it difficult to monitor land cover change and the drivers of change. For instance, the most recent official map of land cover at the national level uses Landsat images ranging from 1999 to 2011 [42]. This and other previous studies have attempted to derive land cover maps based on the use of a limited number of images in areas where cloud cover is frequent. In this study, we overcame this limitation by generating a cloud-free land cover map for the Pacific region of Colombia using all optical data available for the years 2014-2018. Then, by harmonizing this data with other datasets (GFC, UNODC), we were able to quantify and provide understanding of drivers of forest change. However, caution must be taken as the uncertainty of the current land cover increases as the number of years included in the composite increases.
The sample design used for the collection of reference data for land cover classification and validation can greatly influence the robustness of accuracy estimates. In order to minimize any bias in accuracy estimates associated with the collection of training and validation data, both the location of the reference points and their subsetting into training and validation were performed through random sampling. However, we did not test for spatial autocorrelation between sampling locations, which could further influence accuracy estimates [59]. Although autocorrelation between calibration and validation data might occur, we do not see any reason to believe that autocorrelation might substantially influence accuracy estimates or the overall trends. Forest cover is the most abundant class and is the one with the highest likelihood to express strong spatial autocorrelation. We calculated the distance between training and validation points of the most important class, the forest class; from 1969 training points 102 are within 100 m of validation points, 542 are within 500 m of

Discussion
Official land cover maps are not regularly updated in Colombia, making it difficult to monitor land cover change and the drivers of change. For instance, the most recent official map of land cover at the national level uses Landsat images ranging from 1999 to 2011 [42]. This and other previous studies have attempted to derive land cover maps based on the use of a limited number of images in areas where cloud cover is frequent. In this study, we overcame this limitation by generating a cloud-free land cover map for the Pacific region of Colombia using all optical data available for the years 2014-2018. Then, by harmonizing this data with other datasets (GFC, UNODC), we were able to quantify and provide understanding of drivers of forest change. However, caution must be taken as the uncertainty of the current land cover increases as the number of years included in the composite increases.
The sample design used for the collection of reference data for land cover classification and validation can greatly influence the robustness of accuracy estimates. In order to minimize any bias in accuracy estimates associated with the collection of training and validation data, both the location of the reference points and their subsetting into training and validation were performed through random sampling. However, we did not test for spatial autocorrelation between sampling locations, which could further influence accuracy estimates [59]. Although autocorrelation between calibration and validation data might occur, we do not see any reason to believe that autocorrelation might substantially Remote Sens. 2020, 12, 1235 11 of 16 influence accuracy estimates or the overall trends. Forest cover is the most abundant class and is the one with the highest likelihood to express strong spatial autocorrelation. We calculated the distance between training and validation points of the most important class, the forest class; from 1969 training points 102 are within 100 m of validation points, 542 are within 500 m of validation points and 820 are within 1000 m of validation points. The rest, 1149 training points have a distance above 1 km from the validation points. In absence of a rigorous assessment, we advise caution about the potential influence of autocorrelation in some of the results obtained here.
Forests and wetland forests, corresponding to 67,473 km 2 (68%) and 7282 km 2 (7%) of the total area, indicate that large extents of forests remain intact in this biodiversity hotspot. But intensive land use, where native vegetation has been completely removed are present to the north as cattle grazing (1% of the area) and to the south as industrial cropland (1% of the area). The number of deforested pixels that belong to the current grassland class were two fold those from croplands during the study period, but after 2017 the number of deforested pixels associated to the cropland class was larger than the deforested pixels associated to grasslands, indicating a possible change in the trend of these deforestation drivers. Annual deforestation values are especially high from 2006 to 2011, a period in which environmental protection agencies of the Pacific reported that 58.4% of the national timber was produced from natural stands [60]. It is also important to note the large effect that deforestation had during 2016-2017 on wetlands (Table 5), especially in 2016 when 123 km 2 were deforested. Meyer et al. [4] estimated that from 97,024 km 2 of "terra firme" forests, 71,715 km 2 are intact forests, 11,475 km 2 have light degradation, 9234 km 2 have moderate degradation and 4599 km 2 have severe degradation.
A total of 2324 km 2 were identified as forest loss in the period 2001-2018, where the shrubland class was associated with the highest forest loss, followed by wetlands and grasslands. However, one of the limitations of this study is the inability to track changes of land use after coca crop establishment. A major challenge of tracking this driver of forest loss is the poor spectral separability between coca fields and other land cover types such as natural regeneration or agricultural mosaics. Additionally, this crop is highly dynamic as it involves leaf harvesting, replanting, abandonment and migration to forest areas due to eradication policies.
The degradation of wetlands is a major concern given that it is a drug exit corridor [61]; cocaine has been historically trafficked from Colombia to the United States through Central America and Mexico via the Pacific Ocean where rivers are natural corridors facilitating connectivity between cropping areas and distributors.
Based on UNODC reports and GFC forest loss, we estimate that about 25% of the study area is directly affected or under the influence of illicit crops or alluvial mining and that 60% of forest loss was associated with these two drivers. These results suggest that alluvial mining and illicit crops cause more than half of the forest loss during the study period and that land use dedicated to industrial croplands may be gaining importance over cattle grazing.
The large increase in coca crops extent from 2013 to 2018 was accompanied by the largest increase in deforestation between 2016 and 2017 ( Figure 5). However, this region became an important coca producer since 2002. In 2001, UNODC found that the Putumayo and Caquetá departments, which drain to the Amazon River, had the largest amounts of illicit crops in Colombia, but these illicit crops followed a sharp decrease from 2001 to 2012 due to the drug eradication policy (Rincón-Ruiz et al. 2016). Coca crops in the neighboring department of Nariño in the Pacific region, gained importance during 2002 right after 40.000 ha were eradicated from Putumayo and Caquetá [54].
Colombia is one of the main gold producers in Latin America. The National Mining Agency estimated a national production of 61,805 kg during 2016, with 348 km 2 (42%) of alluvial mining located in the collective territories of the Pacific coast [62]. The spatial extent of areas under the influence of alluvial mining and illicit crops are available from UNODC reports. The important extent of deforestation connected to barren land from 2014 to 2017 is likely associated with alluvial gold mining (Table 5). Negative effects of mining operations include water ponds and barren lands resulting after vegetation removal. A major military and police operation against illegal mining in the department of Chocó during 2018 [63] coincided with a drastic reduction in deforested pixels associated with this driver (Figure 4). There is a large presence of alluvial mining in collective Afro-Colombian territories, especially along the Quito River. Radar data allowed identification of areas affected by alluvial mining while optical data allowed identification of barren lands and ponds left after gold mining (Figure 3). Illicit crops and alluvial mining are dispersed over the study area, however, concentration of coca crops is found in Cauca and Nariño, while alluvial mining is concentrated mostly in the Chocó department. Collective territories of El Cantón de San Pablo and Nóvita were the most affected, with degraded areas of 44 and 48 km 2 respectively [62], which is aggravated by the use of mercury in artisanal and small scale mining [64,65]. Vélez et al. [66] found that land titling programs were an important strategy to reduce deforestation rates, however, Afro-Colombian collective territories are heavily affected by this land use, while indigenous reserves are almost not affected [62]. See Appendix A to access available layers of the project web server.
The spatially continuous land cover map, with its assessment of accuracy, enable us to determine that approximately 75% of the study area is covered by forests. Despite the persistent cloud cover of humid/pluvial forests, we were able to determine land cover classes and associate each class to forest loss from 2001 to 2018, improving our understanding of the fate of forest loss pixels after disturbance, and forecasting the drivers of further deforestation. We also found that more than half of forest loss is due to the combined effect of alluvial mining and illicit crops, providing more information to decision makers and conservationists about the magnitude of deforestation associated with the main drivers of forest loss.

Conclusions
Remote sensing is the only alternative enabling the identification of deforestation and land cover conversion in remote areas, where access for field interviews or data collection is constrained by social conflicts and illegal activities. In addition to commercial logging and the increase of the agricultural frontier, this area has been severely affected by alluvial mining and illicit crops. All these drivers combined led to different levels of degradation in the Chocó region. The large extent of shrublands is a result of misuse of this exceptional ecosystem, but it also reflects the ability of this ecosystem to recover from disturbance. Additionally, the shrubland class offers an opportunity for restoration or for the implementation of activities that can reduce further pressure over forests.
In this study, we demonstrate that the Hansen deforestation map can be combined with other data sources to identify drivers of forest loss beyond the generic "increase of agricultural frontier". There was cumulative deforestation of 2324 km 2 from 2001 to 2018, with a trend toward faster deforestation (positive slope). The condition of deforested pixels from 2001 to 2018 based on the current land use increases our understanding of the fate of forests following disturbance. Determining the location and condition of deforested pixels is important in terms of conservation policies. Current conditions of deforested pixels have implications in terms of the degree of degradation of the land and the ecological strategies for restoration. In addition, it affects the cost of the land, from the low cost of degraded lands due to overgrazing to the high cost of lands dedicated to industrial crops.
Our results are valuable to understand the consequences of policy and trade on forest conversion to other land uses including the use of natural resources, restoration, and the price of commodities, such as the demand for cocaine or gold. This study can inform similar analysis in other countries with areas affected by illicit coca crops and mining such as Peru, Bolivia, and more recently Ecuador.