Analysis and Quantification of the Distribution of Marabou (Dichrostachys cinerea (L.) Wight & Arn.) in Valle de los Ingenios, Cuba: A Remote Sensing Approach

Abstract

Cuba is struggling with a growing environmental problem: the uncontrolled spread of the allochthonous weed species marabou (Dichrostachys cinerea) throughout the country. Over the last 70 years, marabou has become a formidable invasive species that poses a threat to Cuban biodiversity and agricultural productivity. In this paper, we present a free and affordable method for regularly mapping the spatial distribution of the marabou based on the Google Earth Engine platform and ecological surveys. To test its accuracy, we develop an 18-year remote sensing analysis (2000–2018) of marabou dynamics using the Valle de los Ingenios, a Cuban UNESCO World Heritage Site, as an experimental model. Our spatial analysis reveals clear patterns of marabou distribution and highlights areas of concentrated growth. Temporal trends illustrate the aggressive nature of the species, identifying periods of expansion and decline. In addition, our system is able to detect specific, large-scale human interventions against the marabou plague in the area. The results highlight the urgent need for remedial strategies to maintain the fragile ecological balance in the region.

1. Introduction

The global landscape is constantly being shaped by the interplay of various ecological forces, with the invasion of alien species standing out as a powerful agent of change. In this intricate ecological web, allochthonous invasive weed species have proven to be particularly resilient actors that significantly impact sensitive ecosystems [1]. According to the European Commission’s definition of invasive weeds from Regulation (EU) No. 1143/2014, invasive plant species are those that, once introduced or spread, threaten or negatively impact biodiversity and certain ecosystem services. In recent decades, invasive weeds have severely threatened local biodiversity, ecosystem services, environmental quality and human health [2]. The United Nations (UN) Intergovernmental Platform on Biodiversity and Ecosystem Services (IPBES) estimates that about one-fifth of the Earth’s surface, including global biodiversity hotspots, is threatened by biotic invaders. This means that 60% of global species extinctions have been caused by invasive alien species alone or in combination with other causes, at an annual cost of USD 423 billion [3]. In these terms, the total reported cost of invasive species in Central and South America between 1975 and 2020 was USD 102.5 billion, with an average cost of USD 2.2 billion per year [4].

Several hypotheses have been proposed to explain the invasion of invasive weeds into new habitats, such as the release of enemies, new weapons and empty niches. However, not every hypothesis is sufficient to explain the invasion of these species [5]. The enemy release hypothesis states that some invasive weeds are more successful in new habitats when they are kept away from natural enemies such as pathogens and herbivores found in their original habitats, which is the case for marabou [6]. Furthermore, the functioning of ecosystems on islands such as Cuba is more affected by invasive weeds than on the mainland [7].

Various invasive plants introduced for the benefit of humans are known to cause ecological and economic damage as they compete for nutrients, light and water [8] and cause significant losses in agriculture and livestock as they increase the costs of their effective management [9]. Their impact is not only economic, as many of these weeds are known to impair ecosystem services such as aesthetics, recreation, culture and regulation [10]. Invasive weeds also affect regulating ecosystem services such as hazard mitigation (e.g., landslides), water purification, pest control, pollination and climate change, which are inextricably linked to agriculture and forestry [10,11].

Eradicating an invasive weed species is a complex process that usually requires extensive resources and long-term monitoring of invasive plant dynamics. The eradication of weeds is often limited at the local level due to various constraints, such as the removal of weeds that spread over large areas. In this context, the regeneration capacity of important perennial weeds such as marabou is also a limiting factor, requiring additional and repeated treatments [12]. Given the negative impact of many invasive weeds, there is an urgent need to prioritize and develop cost-effective and ecologically sound strategies to control them [13].

Marabou (Dichrostachys cinerea (L.) Wight & Arn.), a shrub or small tree with twisted trunks, grayish bark and smooth spines on the lateral branches, reaches a height of 4 to 5 m, even up to 10 m on suitable soils. It was first introduced to Cuba a century and a half ago. However, it was only in the second half of the 20th century that this seemingly harmless plant developed into a major ecological problem. Its high stress tolerance combined with its numerous thorns and resistance to cutting and burning contributed to its uncontrolled spread. Other reasons for the uncontrolled spread of the marabou plant are the dispersal of its seeds by livestock and the high persistence of these seeds in the soil. In addition, the plant reproduces by root buds and has the remarkable ability to produce shoots even when it is not completely eradicated. The collective impact of these traits is the abandonment of land where marabou occurs, as farmers struggle with the physical and economic challenges associated with eradicating the plant. As a result, this weed has spread massively throughout Cuba, particularly in areas already contaminated by overexploitation and abandoned farmland [14].

The Valle de los Ingenios, a UNESCO World Heritage Site, is witness to the ecological consequences of the uncontrolled spread of the marabou plant. The Valle de los Ingenios is a series of three interconnected valleys located about 12 km outside the city of Trinidad and covering an area of 270 km². These valleys—San Luis, Santa Rosa and Meyer—were a center of sugar production from the late 18th to the late 19th century. At the height of the Cuban sugar industry, more than fifty sugar cane mills were in operation in the three valleys, employing over 30,000 slaves. When the valley was declared a UNESCO World Heritage Site in 1988, it was home to dilapidated sugar mills, summer houses, slave barracks and other facilities associated with sugar cane cultivation, as well as active sugar cane plantations. However, the unstoppable advance of the invasive weed marabou has left an indelible mark on this historically and culturally significant landscape [15].

Given the enormous challenge of field assessment, as it is almost impossible to reach the marabou thicket amidst the ubiquitous bush with its many long, hard thorns, remote sensing techniques have proven to be important tools. These techniques play a crucial role in identifying and classifying the weed in various regions of Cuba, including the Valle de los Ingenios [16,17]. Remote sensing, a pivotal component of our study, leverages advanced aerospace technology and geospatial data analysis to monitor vegetation cover through satellite imagery [18,19]. Image processing techniques, such as digital classification, play a critical role in extracting information from satellite images based on digital layers (DL) or reflectance. This process involves grouping pixels based on their values and associating them with specific land cover types to generate detailed maps and directories [20]. Digital image classification, a form of pattern recognition, uses surface features to recognize patterns associated with the location of each pixel [21,22]. Two different classification methods, supervised and unsupervised, depend on how the training statistics are obtained [21]. In the supervised classification approach used here, the process is driven by terrain-specific prior knowledge. Users delineate representative areas for each category of interest, and the pixels within these areas are used to assign the remaining pixels to categories based on the similarity of their reflectance or DL values. The input data can come from field studies, photo interpretations, articles or maps of the area of interest [21]. As described by Hasmadi et al. [23], the process of digital image classification involves three phases: (i) the training phase, in which the categories of interest are digitally defined; (ii) the mapping phase, in which the image pixels are grouped into categories; and (iii) the testing phase, in which the results are verified. In these terms, remote sensing, especially when integrated across different spatial and temporal scales, has the potential to shed light on the dynamics of invasions and reveal their drivers, which can improve the efficiency of control measures [24].

These remote monitoring techniques are of particular interest for agricultural cooperation/development projects in developing countries, such as Cuba. In these countries, access to the desired areas is not always possible or granted; therefore, low-cost, regular remote monitoring would be an invaluable tool for both international cooperation agencies and local authorities [25]. Unfortunately, neither access to regular, expensive, high-resolution satellite imagery required for such consistent analyses nor expensive GIS software licenses are always available or affordable for local authorities, thus hampering the desired technology transfer between international agencies and local institutions. As a starting point for various scientific and agronomic analyses, the quantitative description of spatial patterns of plant species [26] has improved and accelerated in recent years with the development of new geographic information technologies, such as remote sensing [27,28]. These technologies, using medium- and high-resolution satellite imagery, have already shown the ability to map global forest associations or similar detailed floristic contexts [29,30]. The development of an analytical tool capable of remotely quantifying the degree of marabou infestation in large areas, using free cloud-based geospatial analysis platforms and satellite imagery, would be of great interest for assessing the problem of marabou infestation development in Cuba. To achieve this goal, the Google Earth Engine (GEE) is used in this study. GEE is a cloud-based platform that serves as a repository for various datasets from different satellites and platforms. GEE offers the advantage that no satellite images need to be downloaded, and the algorithms can be applied directly to the products stored on the servers. This not only increases processing speed, but also facilitates the development of applications with scalable analysis capabilities in both space and time. In addition, GEE is a free service that can be run on not very demanding computers with a medium-speed Internet connection.

The analysis of satellite imagery involves a multifaceted approach aimed at assessing the distribution patterns of marabou populations over time. Utilizing field data collected in 2014, a classifier was developed using the Google Earth Engine platform. This classifier forms the cornerstone for evaluating marabou populations in subsequent years. The methodology encompasses the creation of composite images representative of each year, synthesizing mean values from valid pixels across the temporal spectrum. The integration of Landsat imagery, corrected for surface reflectance using the USGS Level 2-A product, coupled with cloud and cloud shadow detection using Fmask, ensures the reliability of the analysis. Furthermore, the study introduces a temporal filter to discern between perennial vegetation, such as marabou, and seasonal vegetation, thereby offering a nuanced understanding of its distribution patterns over time. Employing machine learning techniques, the classification algorithm, trained on the 2014 field data, optimizes accuracy and reliability.

The analysis of images from the Valle de los Ingenios over a period from 2000 to 2018 offers a unique opportunity to evaluate the effectiveness of this new tool. This region has witnessed significant human-induced events linked to the marabou infestation, providing well-documented reference parcels for assessing the tool’s accuracy. Looking forward, the application of this tool holds promise for providing stakeholders, including farmers, agronomists, and local authorities, with a comprehensive understanding of the evolving dynamics of marabou infestations in the Valle de los Ingenios. By examining spatial and temporal patterns, it aims to offer valuable insights into the distribution trends of marabou over time.

The aim of this study is to provide a comprehensive, economically affordable tool to analyze the spatial and temporal distribution of the marabou in this region. Both analytical and spatial approaches will be used to understand the dynamics of its invasion and provide valuable insights to both scientists and land managers concerned with the consequences of invasive species in unique and culturally significant landscapes.

2. Materials and Methods

2.1. Study Area

The experimental model for this study is the Valle de los Ingenios, which is located in the province of Sancti Spíritus, Cuba, and covers an area of approximately 250 km². This region, a UNESCO World Heritage Site since 1988, together with the city of Trinidad, faces significant ecological challenges due to the invasive marabou plant (Dichrostachys cinerea). In addition, the area has been affected by two major human interventions in the last 20 years that have severely impacted the spread of the marabou plant in the area. One intervention was the abandonment of sugarcane processing in the valley, which began in 2002 and ended in 2005 with the closure of Central FNTA, the last sugar factory in the area [31]. This led to both a significant decline in agricultural activities and an increase in abandoned agricultural land. The other milestone was the clearing of large marabou-infested areas in the valley with heavy machinery by the national government in 2014–2015 [32]. These two well-known, time-limited human interventions served as benchmarks to ensure the accuracy of our methods.

The precise delineation of the Valle de los Ingenios was initially challenging due to the lack of accurate cartographic data. Initially, efforts were made on the GDAM servers to delineate the Sancti Spiritus region, but this approach proved untenable as it included coastal and mountainous areas beyond the intended scope [33]. As a result, the University of Sancti Spiritus (Cuba) provided images that were crucial to the study, albeit less detailed and to scale (Figure 1). These images were georeferenced and played a crucial role in the precise delineation of the valley. In the search for an optimal boundary, a systematic assessment of different elevation levels was carried out: 80 m, 100 m and 120 m. The aim was to determine the most suitable altitude for the demarcation at hand. The main objective of this study was to determine the altitude that would allow the most accurate and meaningful delimitation of the identified valley boundaries.

2.2. Methodology

2.2.1. Data Source

Satellite images and image processing: A comprehensive dataset of Landsat 5, Landsat 7 and Landsat 8 satellite images from the period 2000 to 2018 was used for the study. The images with a spatial resolution of 30 m were carefully curated using the Google Earth Engine (GEE) platform. It is important to take into account the operating times of the satellites: Landsat 5 (TM) was active from March 1984 to January 2013, Landsat 7 (ETM+) from April 1999 to the present, and Landsat 8 (OLI) from February 2013 to the present. Despite the advantages offered by the Landsat family, challenges have arisen due to missing data and gaps caused by sensor failures. For example, Landsat 7’s ETM+ sensor has provided images with missing data in the form of intermittent bands since June 2003, making analysis difficult. However, by combining data from different sensors, the temporal resolution could be increased to 7 days, allowing a more detailed understanding of Marabu’s annual cycle. In this study, the optical bands of the electromagnetic spectrum played a crucial role in the analysis. Specifically, bands 1, 2, 3, 4, 5 and 7 of the Thematic Mapper (TM) and Enhanced Thematic Mapper Plus (ETM+) sensors and bands 2, 3, 4, 5, 6 and 7 of the Operational Land Imager (OLI) sensor were used [34]. In addition, three bands relating to the temporal development of the Normalized Vegetation Index (NDVI) will be considered in more detail later.

2.2.2. Data Preprocessing

Satellite imagery and data preprocessing: Landsat imagery was used directly in its surface-reflective form, in particular, the USGS Level 2-A product. This approach was complemented by an atmospheric correction using the MODTRAN model to ensure an accurate radiometric calibration. In addition, a cloud and cloud shadow mask created using Fmask [35] proved helpful in overcoming the challenges associated with heavy cloud cover in regions such as Cuba. Given the dynamic nature of cloud-covered areas, of the 23 possible images resulting from the combination of Landsat observations, only about 6 or 7 valid images could be obtained per year. This complexity was further increased by variations in the valid zones at different times.

2.2.3. Harmonic Curve Development

The harmonic curve method uses a mathematical model to represent the variation in NDVI values over time [36]. The model can be expressed as follows:

NDVI = A + B × cos(2πft + φ)

where NDVI is the normalized difference vegetation index, A is the amplitude, B is the phase, φ is the mean value and f is the frequency. The frequency indicates the number of cycles per year, which varies according to vegetation type and location.

The first step in developing the harmonic curve was to calculate the mean values for each pixel per cycle. This was carried out on a grid with three bands representing the parameters of the harmonic NDVI curve: phase, amplitude and mean. The calculation was performed by interpolation using the Google Earth Engine (GEE) platform, which ensured the accuracy of the temporal dynamics representation. In addition, the NDVI was calculated for each scene, which provides important information about the state and dynamics of the vegetation. To improve data quality and reliability, a NoData filter was applied to systematically remove invalid pixels and artifacts from the dataset.

Harmonic curve integration: In addition, the harmonic curve method involved extracting the phase, amplitude and mean of the annual NDVI curve for each pixel using GEE. These parameters served as concise representations of the NDVI curve features for each cycle. The three parameters of the harmonic NDVI curve were integrated into a grid together with the median values of the individual Landsat bands. This integration resulted in a composite image with a total of nine bands for each cycle, providing a comprehensive dataset for further analysis. This combined approach ensured a robust representation of the evolution of the harmonic curve and combined the strengths of mathematical modeling and remote sensing techniques. The integration of these techniques formed a coherent methodology for obtaining meaningful information on vegetation dynamics over time in the Valle de los Ingenios.

2.2.4. Field Data

The ground truth required to develop the classifier comes from extensive field data collected in 2014, with fieldwork focusing primarily on identifying different vegetation types. Some soil types, such as water and urban areas, were identified directly in the images as there was no possibility of confusion. A total of 29 data parcels were strategically selected to represent different land cover classes, including marabou (Mar), palm (Pal), riparian vegetation (Rib), sugarcane (Caña), king grass (King), cassava (Cass), pasture (Past), banana (Ban), etc., which were meticulously recorded (

Table 1. Sampling plots.

Parcel Number	UTM Coordinate (m) WGS_1984_17_N	Area (m2)	Vegetation
1	X: 604,902 Y: 2,413,342	1800	Sugarcane
2	X: 605,243 Y: 2,413,431	900	Cassava
3	X: 605,129 Y: 2,413,855	1350	Pasture, aroma and marabou
4	X: 604,902 Y: 2,413,342	1800	Sugarcane
5	X: 606,551 Y: 2,413,952	5100	Banana
6	X: 606,354 Y: 2,413,808	3500	Sugarcane
7	X: 615,644 Y: 2,418,961	2100	Sugarcane
8	X: 616,699 Y: 2,419,043	15,600	Banana
9	X: 616,812 Y: 2,418,537	39,000	Banana
10	X: 617,175 Y: 2,417,453	5250	Albizia and marabou
11	X: 615,644 Y: 2,418,961	2100	Guava and mango
12	X: 615,581 Y: 2,417,194	7200	King grass
13	X: 614,521 Y: 2,414,653	10,000	King grass
14	X: 615,151 Y: 2,413,776	62,700	Marabou, albizia, guaban, guarana
15	X: 614,123 Y: 2,416,299	43,700	Marabou, aroma
16	X: 614,911 Y: 2,414,276	22,200	Marabou, aroma
17	X: 615,137 Y: 2,414,448	10,500	Marabou, aroma
18	X: 608,456 Y: 2,415,134	28,000	Marabou, aroma
19	X: 608,207 Y: 2,414,556	24,000	Eucalyptus
20	X: 611,645 Y: 2,415,412	85,800	Marabou, albizia, aroma
21	X: 610,162 Y: 2,415,228	18,000	Marabou, guin de bandera
22	X: 607,502 Y: 2,414,784	36,975	Marabou
23	X: 607,502 Y: 2,413,310	24,000	Marabou and riparian species
24	X: 614,138 Y: 2,414,879	13,500	Marabou, guacima
25	X: 614,143 Y: 2,414,625	12,600	Marabou, aroma
26	X: 613,737 Y: 2,415,050	56,700	Marabou
27	X: 607,813 Y: 2,415,031	9000	Riparian vegetation with some marabou along the edges
28	X: 607,833 Y: 2,413,688	18,000	Marabou, aroma, guarana, palma
29	X: 606,533 Y: 2,412,695	2750	Palm, almacigo, yagruma, guacima

). This diversity of classes was intended to capture the broad spectrum of land cover in the Valle de los Ingenios. To address concerns about the classification process and the purity of the ground truth plots, detailed information has been included in Appendix A. This section provides a detailed description of the training data and sample plots and explains the process used to derive the random forest classifications. It also includes information on the purity of each plot, addressing questions about the presence of various elements such as roads, agricultural infrastructure, residential buildings, etc. In addition, Appendix A serves to reconcile any inconsistencies in the land cover class names in the different sections of the document to ensure transparency and clarity in our classification methodology. A deliberate choice was made to abstain from creating a land-use map, aligning with the study’s core objective of just monitoring changes in marabou distribution within the Valle de los Ingenios. For the training and validation of the model, high-resolution satellite imagery was used to extract multiple pixels within each of the 29 plots. The resulting confusion matrix, comprising 136 pixels, reflects the classification of each pixel within these plots using the trained random forest (RF) model. This approach ensures a comprehensive and representative assessment of the model’s performance.

2.2.5. Classification Algorithm

In the context of this study, the classification criterion employed for the distribution of marabou is straightforward: every pixel that contains any amount of marabou is unequivocally classified as marabou. In this study, the image classification methods in the Google Earth Engine (GEE) were examined in detail. The methodology included a thorough evaluation of four different supervised classification algorithms within GEE: CART, RandomForest, NaiveBayes and SVM. The process involved the creation of confusion matrices for each algorithm with the aim of identifying the most appropriate classifier that could explain the predominant variance in the dataset.

The parameterization of the random forest algorithm (RF) included the selection of the following input variables from the Landsat 8 images: B2 (blue band), B3 (green band), B4 (red band), B5 (near-infrared band), B6 (shortwave infrared band 1), B7 (shortwave infrared band 2), the phase band, the magnitude band and the NDVI (calculated from the red band and the near-infrared band). The number of trees in the RF classifier was set to 100, which was considered an appropriate value for this dataset. However, it is important to note that the optimal number of trees may vary depending on the specific characteristics of the data. To ensure the generalizability of the classifier, it was evaluated on a separate dataset consisting of 30% of the images in the sample region. For transparency and reproducibility, the full programming code for the RandomForest algorithm in Google Earth Engine (GEE) is included in Appendix B. Researchers and practitioners can refer to this appendix to gain a detailed understanding of the implementation and configuration of the algorithm.

2.2.6. Workflow Overview

The workflow overview is summarized in Figure 2.

2.2.7. Temporal Considerations

A classifier developed from the 2014 field data was systematically applied to each annual cycle. For periods prior to 2014, semi-annual cycles were considered due to the limited Landsat 5 data. From 2014 onwards, Landsat 7 and Landsat 8 images acquired every 8 days enabled annual classification. The application of the classifier resulted in maps that illustrate the temporal evolution of marabou vegetation in the Valle de los Ingenios for each cycle and provide valuable insights into vegetation dynamics over time.

2.2.8. Code Accessibility

To ensure transparency and reproducibility, the entire code that implements the classification and analysis processes is openly accessible via the GEE code editor at code.earthengine.google.com (accessed 1 August 2023). This accessibility, which is described in detail in Appendix B, facilitates verification and replication of the methodology, thereby increasing the robustness and verifiability of the results. Researchers and practitioners interested in exploring or validating the approach can easily access the code to gain a comprehensive understanding and application.

3. Results

3.1. Valley Mapping at an Elevation Threshold of 100 m

The integration of digital terrain models [37] (DTM) with satellite imagery was crucial for the delineation of the Valle de los Ingenios. By delineating the valley, we aim to establish clear spatial boundaries for our study area, which is crucial for accurate analysis and interpretation of the data. After testing different elevation levels, the application of an elevation threshold of 100 m made it possible to delineate the northern, eastern and western extent of the valley. The overlay of these boundaries with the images provided by the local authorities resulted in a coherent representation of the topography of the valley (Figure 3).

However, the southern boundary proved to be a difficult task as elevation 100 alone could not fully capture the intricate features of the valley. To solve this problem, various mountain ranges bounding the southern edge of the valley were considered in conjunction with high-resolution image interpretation. Thank to this comprehensive strategy, the southern boundaries could be precisely determined, ensuring an accurate representation of the Valle de los Ingenios in its entirety.

The resulting map, derived from the integration of elevation data and satellite imagery, provided a visually compelling representation of the geographic contours of the valley. Our digital delineation yielded a total area of 22,919.12 ha, which is close to the 250 km² normally attributed to the valley. This elevation-based mapping approach served as a fundamental step in understanding the spatial characteristics of the Valle de los Ingenios and formed the basis for subsequent analysis and interpretation.

3.2. Harmonic Curves

Harmonic curves were used to represent the temporal variability of the NDVI values for each pixel of the study area. The curves were calculated using the integrated tools of Google Earth Engine (GEE) and are based on data from Landsat 5, Landsat 7 and Landsat 8 satellite imagery from 2000 to 2018 (Figure 4).

The results of the harmonic curve analysis provided insight into the temporal dynamics of marabou vegetation in the study area and showed that marabou vegetation in the study area exhibits a seasonal cycle with peaks in NDVI during the local lighter, wet season (June to October) and minima during the darker, dry season (November to May) (Figure 5). This is consistent with a phase of 0.75 indicating that the peak in NVDI values occurs in the middle of the lighter, wetter season (Figure 5). These harmonic temporal NDVI changes are common in perennial species and may be related to changes in the photosynthetic activity of marabou due to the alternation between wet and sunny or dry and dark conditions [38]. The mean NDVI value for marabou vegetation in the study area was 0.75, indicating that marabou is a relatively productive vegetation species. The mean NDVI value for marabou vegetation in the study area was calculated to be 0.75, with an amplitude of 0.48. This suggests that NDVI values exhibit variability of up to 40% over the course of a year, which is consistent with NDVI variability in tropical vegetation, where sunlight and water levels can fluctuate over time but are not limiting factors [39].

3.3. Confusion Matrix

Based on field data collected in 2014, the random forest algorithm proved to be the best choice, exhibiting good overall accuracy. With an overall accuracy of 0.92, Random forest categorized more than 92% of the data correctly and showed the best accuracy in assigning labels to different land cover classes. Furthermore, the kappa coefficient yielded a value of 0.89. This robust kappa value shows that the random forest algorithm was able to recognize patterns within the dataset, even when randomness was taken into account (

Table 2. Confusion matrix from the 2014 random forest classification.

Class	Mar	Pal	Rib	Cane	King	Urb	Water	Mang	Sum
Mar	56	0	0	0	1	0	0	0	57
Pal	0	9	0	0	1	0	0	0	10
Rib	0	0	2	0	0	0	0	1	3
Cane	0	0	1	17	0	0	0	1	19
King	1	1	0	0	3	0	0	0	5
Urb	0	0	0	0	0	34	0	2	36
Water	0	0	0	0	0	0	4	1	5
Mang	0	0	0	0	0	0	0	1	1
Sum	57	10	3	17	5	34	4	6	136

OA: 0.9265. Kappa: 0.8996.

).

3.4. Coverage Maps and Data

The analysis of the development of the marabou infestation in the Valle de los Ingenios revealed a complex pattern of spatial distribution and temporal development. Originally, the marabou area was relatively small, covering 7908.1 hectares in 2002, until it declined to 5067.5 hectares in 2004. During this period, the situation remained stable, and the infestation was mainly concentrated in the northern, eastern and western borders of the valley, bordering the Sierra del Escambray, a mountainous region with limited agricultural activity and difficult access (Figure 6a,

Table 3. Evolution of marabou infestation in the Valle de los Ingenios area between 2002 and 2021 year.

Year	Hectares
2002	7908.1
2004	5067.5
2006	14,376.8
2008	14,132.2
2010	12,326.2
2012	12,945.2
2014	12,670.7
2015	12,796.1
2016	10,330.2
2017	8772.5
2018	9329.7
2019	12,675.0
2020	15,657.0
2021	15,640.0

).

Between 2004 and 2005, coinciding with the closure of the last sugar factory in the region (Central FNTA or Central Trinidad), the area affected by marabou increased dramatically, tripling from 5067 hectares in 2004 to an astonishing 14,376 hectares in 2006. A comparison of the distribution maps for these years shows that the marabou’s main areas of distribution are along the Sierra del Escambray, which has served as a weed reservoir since 2002, and along the Carretera 12 road (which connects Trinidad to the town of Sancti Spiritus), which has proven to be an important distribution route for the weed (Figure 6b, Table 3).

From 2006 to 2015, the area affected by marabou remained relatively stable, fluctuating between 14,376 hectares in 2006 and 12,326 hectares in 2010 (Appendix C, Table 3). In 2015, coinciding with the Cuban government’s institutional campaign against marabou on the occasion of the 500th anniversary of the city of Trinidad, a significant decrease in the area infested by marabou was observed. The affected area decreased from 12,796 hectares in 2015 to 10,330 hectares in 2016 and further to 8772 hectares in 2017. The cleared areas were mainly located along the main roads, which was due to the accessibility of these areas for heavy machinery and the municipality’s increased focus on clearing areas of high tourist importance (Figure 6c). As more areas were cleared, farmers gained easier access to previously inaccessible areas, which may explain the decrease in marabou-affected areas in 2017 after the government initiative was completed (Appendix C, Table 3).

In the last third of the study period (2018–2021), a progressive increase in the marabou-affected area was observed, rising from 9329 hectares in 2018 to over 15,000 hectares in 2020 and 2021, marking the highest values in the analyzed time series (Figure 6d Appendix C, Table 3). This resurgence can possibly be attributed to the clearing of previously abandoned sugarcane fields that had acted as barriers to the spread of marabou. The removal of these barriers facilitated the spread of the marabou plant into new areas, leading to the observed increase in infestation.

4. Discussion

The rapid spread of the marabou plant, an invasive plant species, poses a significant ecological and economic threat to the Valle de los Ingenios, a UNESCO World Heritage Site in Cuba. In our study, Google Earth Engine (GEE), a cloud-based platform for geospatial analysis, was used to comprehensively assess the spatial and temporal distribution of marabou infestation in the valley. Our results first highlight the importance of using cost-effective and user-friendly tools such as Google Earth Engine to assess marabou infestation in developing countries. GEE’s ability to accurately locate and monitor marabou infestations over time, as well as its ability to detect the impact of human activities on the environment, make it a valuable tool for environmental management and decision making [40,41]. Furthermore, our results are consistent with previous findings that emphasize the importance of time series analysis, particularly the NDVI, for understanding the successional trajectories of both tropical forests and coastal areas [39,42]. The intricate relationship between forest and coastal structure and NDVI can be applied to our context and helps characterize marabou infestation stages in mixed environments such as Valle de los Ingenios.

Our analysis of nearly two decades of satellite imagery unveiled a nuanced pattern of marabou expansion and decline, shaped by diverse factors. Initially, the marabou population remained modest; however, following the closure of the area’s last sugar factory in 2005, it underwent a substantial surge. This surge is attributed to the abandonment of agricultural land consequent to the sugar industry’s decline, fostering ideal conditions for marabou dispersal. The observed dispersal pattern underscores a multifaceted interplay between ecological and socio-economic factors, with abandoned agricultural land facilitating marabou settlement and dispersal. This colonization dynamic adheres to established patterns, including initial dispersal in natural weed reservoirs like hard-to-reach mountainous areas on the valley’s periphery, and the utilization of transportation infrastructure as dispersal routes, common strategies for both allochthonous and autochthonous weed dispersal [43].

Despite the general stability observed the marabou infestation during the 2006–2015 period, the local community actively sought to control the infestation, focusing on accessible areas of high tourist value, such as the main roads, which previously favored the spread of marabou.

The subsequent government-led campaign against marabou in 2015, which coincided with the 500th anniversary of the city of Trinidad, led to a significant reduction in the area infested with marabou. Since then, however, the marabou population has rebounded, possibly due to the clearing (but not subsequent utilization) of previously abandoned sugarcane fields, which previously acted as a barrier to the spread of marabou. As with many other invasive weeds, allochthonous marabou is more competitive than native species and rapidly colonizes cleared areas if they are not restored or used for agriculture. This aggressive colonization dynamic has already been observed in invasive weed species [43,44,45] and poses a real threat to both agriculture and biodiversity. In this respect, our results are consistent with other studies on invasive weed species which are based on phenology-based mapping, such as in Labonté et al. [46]; use satellite images to determine the influence of disturbance and management practices on the spread of invasive plants [47]; or demonstrate the importance of time-series analysis in invasive species mapping [48].

Acknowledging the limitations of our study, particularly in the validation process, is essential. Validation, conducted solely for the year 2014, potentially introduces discrepancies due to differences between satellite platforms. The projection of the 2014 model without further validation presents a limitation; however, the relatively high accuracy achieved using only NDVI suggests that integrating additional complementary parameters could enhance future analyses comprehensively. Variations in satellite platforms may introduce discrepancies that necessitate consideration when interpreting the results. Addressing these limitations is pivotal for improving the accuracy and robustness of invasive species analyses, rendering them more applicable for effective environmental management and decision-making. Nonetheless, it is noteworthy that marabou, being a woody, evergreen plant, exhibits less spectral change compared to annual or deciduous woody plants, suggesting that NDVI may serve as a more reliable indicator of marabou presence than for other types of invasive species.

Our results emphasize the importance of the future evaluation of different marabou control methods, especially given the observed resurgence of marabou after clearing activities. Tailored strategies may be needed that are adapted to Cuban specificity and take into account the specific characteristics of different areas. For example, the promotion of agricultural activity and the reforestation of native forests near the Escambray Mountains and the enhancement of the tourist landscape along the transportation routes in the valley itself could be promising approaches. Our results emphasize the need for a comprehensive approach to marabou management that addresses both the ecological and socio-economic drivers of marabou expansion. This includes measures to restore abandoned agricultural land, the promotion of sustainable agricultural practices and the involvement of the local population in management measures.