iForest - Biogeosciences and Forestry

iForest - Biogeosciences and Forestry

Afforestation monitoring through automatic analysis of 36-years Landsat Best Available Composites

iForest - Biogeosciences and Forestry, Volume 15, Issue 4, Pages 220-228 (2022)
doi: https://doi.org/10.3832/ifor4043-015
Published: Jul 12, 2022 - Copyright © 2022 SISEF

Research Articles

The study of afforestation is crucial to monitor land transformations and represents a central topic in sustainable development procedures, in terms of climate change, ecosystem services monitoring, and planning policies activities. Although surveying afforestation is important, the assessment of the growing forests is difficult, since land cover has different durations depending on the species. In this context, remote sensing can be a valid instrument to evaluate the afforestation process. Nevertheless, while a vast literature on forest disturbance exists, only a few studies focus on afforestation and almost none directly exploits remote sensing data. This study aims to automatically classify non-forest, afforestation, and forest areas using remote sensing data. To this purpose, we constructed a reference dataset of 61 polygons that suffered a change from non-forest to forest in the period 1988-2020. The reference data were constructed with the Land Use Inventory of Italy and through photointerpretation of orthophotos (1988-2012, spatial resolution 50 × 50 cm) and very high-resolution images (2012-2020, spatial resolution 30 × 30 cm). Using Landsat Best Available Pixel composites time-series (1984-2020) we calculated 52 temporal predictors: four temporal metrics (median, standard deviation, Pearson’s correlation coefficient R, and slope) calculated for 13 different bands (the six Landsat spectral bands, three Spectral Vegetation Indices, and four Tasseled Cap Indices). To verify the possibility of distinguishing afforestation from non-forest and forest, given the differences between them can be minimal, we tested four different models aiming at classifying the following categories: (i) non-forest/afforestation, (ii) afforestation/forest, (iii) non-forest/forest and (iv) non-forest/afforestation/forest. Temporal predictors were used with random forest which was calibrated using random search, validated using k-fold Cross-Validation Overall Accuracy (OAcv), and further using out-of-bag independent data (OAoob). Results illustrate that the distinction of afforestation/forest reaches the largest OAcv (87%), followed by non-forest/forest (83%), non-forest/afforestation (75%) and non-forest/afforestation/forest (72%). The different OA values confirm that the difference in photosynthetic activity between forest and afforestation can be analysed through remote sensing to distinguish them. Although remote sensing data are currently not exploited to monitor afforestation areas our results suggest it may be a valid support for country-level monitoring and reporting.

Afforestation, Remote Sensing, Land Cover Monitoring, Random Forest


According to Food and Agriculture Organization (FAO) definition, the forest is any “territory with arboreal coverage greater than 10% compared to an extension greater than 0.5 ha, where the trees reach a minimum height of 5 m when ripe”, including also: “[…] windbreak barriers and wooded strips of a width exceeding 20m, rubber tree plantations, and cork trees” ([10]). Thirty-one percent of the World’s total surface is occupied by forests (4.06 billion hectares) and Europe occupies a prominent position in this context since European forests account for 1,017,461 ha ([11]). Forests play a strategic role in terms of landscape, history, culture, and economy ([40]). They contribute significantly to human well-being, providing a multitude of ecosystem services among which wood production, carbon sink, food production, soil protection, climate control, biodiversity conservation, environmental decontamination, and hydrological cycle regulation ([47]).

The Italian forest heritage is one of the most important in Europe - 11.8 million ha or 40% of the Italian land ([34]) - but it is exposed to many threats, in the majority of cases directly or indirectly related to anthropic activities. Indeed, great pressure is currently given by climate change, which undoubtedly led to an increase in the occurrence of extreme events such as storms, fires, and droughts ([13]). Therefore, it results of primary importance to efficiently monitor forest ecosystems and their status, as well as to detect change processes and provide new tools for their conservation.

Forest change processes are territorial transformations of primary importance, and have a central role in the definition of policies aimed to preserve forest resources and its ecosystem services provision. International organizations and other institutions at different levels have undertaken forest changes issues, giving definitions to frame this phenomenon.

The United Nations Food and Agriculture Organization Forestry Department established a series of definitions related to changes between forest and other land use classes, including the loss of forest cover. According to these definitions, it is possible to distinguish two different processes of tree cover growth: (i) afforestation, i.e., the conversion from other land uses into the forest, or the increase of the canopy cover to above 10%; and (ii) reforestation, i.e., the re-establishment of forest formations after less than 10 years with less than 10% canopy cover due to human-induced or natural perturbations. In this sense, reforestation is presented as a short-term process that, if longer than 10 years, falls within the definition of afforestation ([10]).

Forest changes occur heterogeneously around the globe, mainly because of different socio-economic backgrounds. While in Africa and South America and Eastern Europe deforestation trends up, in the European Mediterranean areas afforestation, generally resulting from agricultural to forest land-use changes, represents broadly the dominant process ([33], [21], [37]). Regarding the Italian territory, historical cartographic and inventory data confirms this assumption, showing that forest surface has been constantly increasing through the years, with an average annual growth rate of 0.3% ([2]).

For the above reason, afforestation represents a process of increasing importance, which certainly needs to be characterized to make the most of its potential and to correctly drive its development in line with forest sustainable management principles. As on one hand, afforestation certainly involves all the benefits connected to forest ecosystem services, on the other hand, if not well managed, it can cause land fragmentation and other related issues ([50], [8]). In other words, most of the issues and challenges associated with afforestation depend on its location ([33]).

Since the afforestation process has a significant effect on the land cover dynamics, it is necessary to monitor its evolution, and to this purpose many different independent data in national and European environments are available, both inventory and cartographic. One of the most important is the Land Use Inventory of Italy (IUTI), an inventory data based on a sampling system and the classification in 6 macro-categories ([32]). Although IUTI is a valuable source of information, it is seldom updated (1990, 2000, 2008, 2013, and 2016) and provides information on the predominant land cover in the area of analysis (in Italy 1,217,032 random points distributed in as many quadrants of 25 ha, with a minimum mapping unit of 5000 m2 and a minimum width of 20 m).

In this context, Remote Sensing (RS) certainly represents a high potential instrument ([29]. [12]), which allows to analyse and classify large areas, and enables to obtain frequently updatable results at different spatial scales ([1]). Remote sensing data analysis can be improved thanks to cloud platforms as Google Earth Engine® (GEE - [19]), which allow managing large amounts of data and information even with small computational capabilities. Google Earth Engine® provides numerous satellite missions as Landsat and the more recent European Sentinel Copernicus program which are by far the most used, being open-source data and covering many years of observations at different spatial resolution ([27]).

Information about land cover characteristics can be evaluated using the different bands of multispectral data, which can be summarised in single values, called Spectral Indices, or, when they are related to vegetated land cover, Spectral Vegetation Indices (SVI - [48]). The most used SVI in literature are the Normalized Difference Vegetation Index (NDVI), the Enhanced Vegetation Index (EVI) and the Normalized Burn Ratio (NBR). Besides the SVI, other indices known as Tasseled Cap (TC) represent a valuable tool for this kind of analysis. These indices resume multispectral data into four dimensions of the land cover considered: Brightness (B), Wetness (W), Greenness (G) and Angle (A). SVI and TC acquire even more value when referred to different periods, allowing to build temporal trends which represent important additional information about the phenology of a given territory ([38]). Hence, to elaborate a land cover classification and to produce forest cartography, an efficient method consists in analysing multitemporal series of SVI or TC in order to attain synthetic though very complete information. In doing so, it is possible to establish classification rules by imposing thresholds to discriminate different land cover classes and recognising the ones of interest ([46], [47], [31], [9]).

Landsat and Sentinel images have different applications related to their intrinsic characteristics. Landsat data are available since the 1970s, thus providing very long time data series at medium spatial resolution (30 m) which are useful to analyse land cover trends and evolution. Instead, Sentinel is available since 2015, providing images with high spatial resolution (10 meters) and short revisiting time (in Italy 2/3 days), which make it a tool of great utility for applications where high precision is needed. Sentinel-1 data and Sentinel-2 images, thanks to their spatial and temporal resolutions, have been recently used in many studies on the analysis of forest disturbances and land cover change detection. In Laurin et al. ([30]), Sentinel-1 data were used to detect the damages caused by the Vaia wind storm in 2018 in Northern Italy, while Sentinel-2 data were used to map clearcuts in Tuscany ([17]), and to estimate areas and relative confidence intervals of forest fires, wind damages and clearcuts in Italy for 2018 (Francini et al., in review).

Other research activities have taken advantage of the high temporal resolution of some satellite images to monitor land cover and land cover changes at the annual or sub-annual level. De Fioravante et al. ([9]) integrated the Sentinel-2 and Sentinel-1 data to assess land cover (2018) and land cover changes (2017-2018) in Italy, while Francini et al. ([15]) used PlanetScope imagery to develop a methodology to detect near-real-time forest disturbances in Tuscany (Italy). While the land cover and its changes monitored by remote sensing is a widely analyzed topic, afforestation assessment remains a theme mainly unrepresented, being treated as a central subject only in a few studies ([41], [44]) in which relatively low accuracy was reached (less than 75%). For example, Huang et al. ([26]) and Qiu et al. ([39]) used MODIS images (500 m spatial resolution) to map forest growth areas in China. However, products with a spatial resolution of 500 meters are not useful for monitoring afforestation in fragmented and morphologically complex areas typical of the Italian territory, where afforestation processes are characterised by patches of limited size.

The present research aims to illustrate a new methodology that exploits the Landsat time series to demonstrate that RS can provide valuable support not only for monitoring stable land cover classes ([9]) and forest disturbances ([18]), but also for predicting afforestation areas.

In the first section of this paper, we describe the training data and Landsat images (1984-2020). Then we explain the use of the input data to assess the temporal trend difference of photosynthetic activity of the classification periods (non-forest, afforestation and forest). Afterwards, we focus on the use of the Random Forest™ (RF) algorithm to automatically distinguish the classes and identify the combination of predictors which allows reaching the largest accuracy.

In the last part of the paper, we discuss the results and how they support the use of remote sensing data in the distinction between non-forest, afforestation and forest.


Study area and reference data

We selected sixty-one areas distributed along the Italian peninsula, especially on the mountains range (Alps and Appennines) where afforestation occurred in the period 1988-2020. Forty-two sites are located in areas with elevation higher than 600 m a.s.l., and twenty-six of them are located in areas above 1000 m a.s.l; the slope of almost 3/4 of the points is higher than 10% (Fig. 1). These areas cover a surface of 269.197 ha, and each polygon is on average 4 ha large. The areas were defined using the IUTI dataset and other images, as explained below.

Fig. 1 - Study area and distribution of reference dataset. On the right, an example of afforestation occurred between 1990-2000 is displayed: in the upper image (1988) there was no forest cover, in the centre (2000) the forest partially cover the area, in the last image (2019) the forest cover is complete.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Land Use Inventory of Italy (IUTI)

Land Use Inventory of Italy IUTI was realized by the Italian Ministry of Environment and Protection of Land and Sea in 1990 (updated in 2000, 2008, 2013, and 2016) in the framework of the Extraordinary Plan of Environmental Remote Sensing. It is composed of 1,217,032 randomly selected points, which have been classified through photointerpretation, considering a minimum mapping unit of 5000 m2 and a minimum width of 20 m. The classification system is based on the Intergovernal Panel on Climate Change guidelines, and it is composed of three hierarchical levels, which aim to identify the land cover categories which are important for the Kyoto Protocol and to integrate the National Inventory on Forest and Carbon Pools (INFC) results on woods and other woody areas categories, based on the FAO definition ([42]). Analysing the evolution of the dataset through time makes it possible to identify those areas where afforestation occurred.

Landsat data and Best Available Pixel Composite

In this study, the Landsat surface reflectance data available (for Landsat mission 5-7-8) in the Google Earth Engine® (GEE) archive were used. We used Landsat missions 5, 7 and 8, because the product of the previous missions had a lower spatial and temporal resolution. Moreover, we analysed data starting from 1985 because they are comparable with the training dataset. Landsat images are composed of six bands: three visible (blue, green, and red), one near-infrared (nir), and two short-wave infrared (swir1 and swir2). All of them are already processed to orthorectify surface and brightness (the thermal infrared) reflectance, atmospherically corrected, and cloud masked. Landsat images were used to create Italian annual composites from 1984 to 2020, from June 1st until August 31st, using the Best Available Pixel (BAP) procedure ([51]).

The BAP aims to fill the final image mosaic with the composite best available pixel surface reflectance value. The selection of the best pixel of the images collection takes into account four different criteria: (i) sensor score, to penalize Landsat 7 images where the Scan Line Corrector malfunction (SLC-off) is present; (ii) target day score, to preferably select the images acquired close to a defined acquisition day (in this case 15th of August); (iii) distance to cloud/cloud shadow score, to decrease the scores of those pixels which are in the proximity the cloud cover; (iv) opacity score (calculated using opacity band produced by LEDAPS - [43]), to prefer the pixel with low opacity. Scores (i) and (ii) are applied to the whole image, while scores (iii) and (iv) are applied to each pixel. BAP composites calculation was performed using the GEE application developed by Francini et al. ([16]) which allows calculating and downloading BAP composites ([24], [25], [51]). For more information on the criteria involved in the BAP pixels selection see the BAP-GEE documentation ([14]). For more info on the BAP algorithm see White et al. ([51]) and Griffiths et al. ([20]).


The methodology presented here employed Landsat images to define the feasibility of the assessment of afforestation areas in the last thirty-six years through remote sensing. The analysis of many years allows verifying the evolution of the classification periods, considering the afforestation, which can have a different duration depending on the species.

Using 61 training areas consisting of polygons defined using IUTI elements, orthophotos and very high resolution images to identify the afforestation, a preliminary analysis of photosynthetic activity in non-forest, afforestation and forest areas was developed. Then RF algorithm, calibrated and validated, was applied to automatically distinguish the non-forest, afforestation and forest land cover (Fig. 2). The analyses were carried out using QGIS, GEE® and R softwares, which are useful in processing spatial and statistical data.

Fig. 2 - Workflow of the applied methodology. The steps to distinguish among non-forest, afforestation and forest areas are briefly described.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Reference data acquisition

The orthophotos used for the photointerpretation process are available as Web Map Services (WMS) of the Italian Ministry of Environment, Land and Sea, and they are referred to 1988-1989, 1994-1998 (black and white) and to 2000, 2006, 2012 (colour), with a spatial resolution of 50 cm. The more recent Very High Resolution images (30 cm spatial resolution, available since early 2000s) are freely available thanks to the web service of QGIS and Google Earth Pro®.

The reference dataset (Fig. 1) was produced through photointepretation of orthophotos and Very High Resolution images, available from 1988 to 2020; in addition to known areas, the IUTI dataset has been used to establish the afforestation areas by considering the points which identify a forest change (from non forest to forest area). Sixty-one training areas have been collected and finally integrated with other information on observation date and polygon dimension.

Orthophotos and IUTI datasets were used to define the temporal range in which the afforestation process occurred. More specifically, two years were identified: (i) the year corresponding to the transition between non-forest and afforestation; and (ii) the year corresponding to the transition between afforestation and forest.

Temporal predictors calculation and analysis

The indices calculation for each composite BAP and the analysis of the resulting temporal series were carried out to calculate the temporal statistics. More specifically, we augmented the six Landsat bands of our BAP composites (see above) by calculating seven additional indices: (i) Normalized Difference Vegetation Index (NDVI); (ii) Normalized Burnt Ratio (NBR); (iii) Enhanced Vegetation Index (EVI); (iv) Tasseled Cap brightness (B); (v) wetness (W); (vi) greenness (G); and (vii) angle (A). As a result of this step, we obtained 13 temporal series of 36 years, hereinafter called Time Series (TS - Fig. 3).

Fig. 3 - Example of NDVI series in an afforestation pixel. For each classification period the trend line was obtained and used to calculate the Pearson’s correlation coefficient (r).

  Enlarge/Shrink   Download   Full Width  Open in Viewer

For each polygon and for each of the non-forest, afforestation, and forest temporal ranges, we calculated from the TS four temporal statistics: (i) the mean and (ii) the standard deviation of the TS; (iii) the slope; and (iv) the Pearson’s correlation coefficient (r) of the regression line obtained by linear interpolation of the TS over time. As a result of this step, a set of 52 temporal predictors (13 TS × 4 temporal statistics) were obtained to distinguish photosynthetic activity between non-forest / afforestation / forest areas and to investigate the differences in photosynthetic trends using density plots.

Random Forest

Random Forest model and hyperparameters

To classify forest, non-forest, and afforestation areas we used a random forest model and the temporal predictors as variables. Random Forest™ (RF) is a decision tree algorithm often applied to forest classification ([23]), where a subset of predictors is used to build a set of regression trees, applying the Out-Of-Bag sample (OOB) procedure ([6]). RF is also useful as it allows to assess the variable importance in the model. RF algorithm was chosen because it achieves large accuracies when compared with other machine learning approaches ([35]) and it is not influenced by overfitting ([4]). Random Forest hyperparameters are model parameters that can be calibrated to obtain the minimum error rate ([30]). In this study, we tuned the following two: max features (the number of features considered) and max depth (the maximum depth of the tree, the longest path between the root node and the leaf node). The number of trees in the algorithm was set to 300 and the minimum number of samples for splitting a node was set to one.

Max depth and max features calibration

Max features and max depth hyperparameters were calibrated using the procedure named random search ([5], [30]). Specifically, we defined a grid of max features (1-52) and max depth (1-40) ranges, from which a subset of 100 combinations were randomly selected. To assess the performance of each tested combination we used the k-fold Cross-Validation (CV) with k = 5.

Cross Validation serves to train and validate the model on different data and it is useful to reduce the overfitting ([52]) and to test the robustness of the model ([30]). Using k-fold CV, the dataset is divided into k groups, which are composed of the same quantity of data: one of these folds is defined as the validation group, while the remaining k-1 groups are considered as training data. The procedure is then repeated k times, each time using a different fold as a validation group. Finally, the validation data is merged and the k-fold CV Overall Accuracy (OAcv) is calculated. As a result of this step, we obtained the max features - max depth combination that allowed to obtain the larger OAcv.

Performance assessment

Although the k-fold CV procedure (see above) allowed to rigorously estimate the accuracy (OAcv), we further assessed the accuracy of the model exploiting the best max features - max depth combination using Out Of Bag (OOB) data, which allows getting an unbiased estimate of the error ([6]). During random forest training, two independent sets are created: (i) the bootstrap sample, i.e., the data chosen to be “in-the-bag” by sampling with replacement; and (ii) the out-of-bag (OOB) sample, i.e., all data not chosen in the sampling process. During model training, this process is repeated and many OOB samples are created. Since they are never used to train the model, they are never-seen-before data and can thus be used to assess the performance of the model. Finally, all OOB samples are aggregated and OOB predictions are compared with reference data by constructing a confusion matrix from which we calculated the OOB. Overall Accuracy (OAoob), omission errors (omissionsoob), and commission errors (commissionsoob). The procedure was repeated four times, to compare the classification models as follows: (i) non-forest/afforestation, (ii) afforestation/forest, (iii) non-forest/forest, and (iv) non-forest/afforestation/forest. For each of them, classification omission and commission error were calculated, to respectively identify the pixel non-classification and the misclassification.

Variables importance calculation

To identify temporal predictors that mostly contributed to increasing the accuracy of the model we analysed the variable importance ranking outputted by random forest. The variable importance was expressed in terms of the Mean Decrease of the Gini index (MDG - [36]). MDG is the measure of the variable contribution to the homogeneity of the random forest trees component. Large MDG values correspond to variables that relevantly contributed to increasing the performance of the model and thus to the accuracy achieved. Conversely, small MDG values correspond to variables that did not contribute to increasing the accuracy of the model.


Classification accuracies

Random forest validation using Out-Of-Bag data (OAoob) was performed on all the temporal predictors of the best iteration of 5-fold CV (Tab. 1). The largest OAoob was reported by the model classifying between afforestation and forest (80%), followed by the classification between non-forest and forest (77%). Non-forest and afforestation and non-forest/afforestation/forest classification model reached an overall accuracy smaller than 70%.

Tab. 1 - Classification models overall accuracy of random forest (OAoob), considering all the temporal predictors.

Classification models max features max depth OAoob
non-forest/afforestation 24 21 0.70
afforestation/forest 7 5 0.80
non-forest/forest 52 30 0.78
non-forest/afforestation/forest 24 28 0.67

  Enlarge/Reduce  Open in Viewer

For each class of the four classification models, omission (omissionsoob) and commission (commissionsoob) have been calculated (Tab. 2). Both omission and commission errors are larger than 10%, with the maximum value (40%) registered by non forest class in the classification between non-forest, afforestation and forest; the smallest error values are registered in the afforestation and forest classification, where omission error of the forest class is 11% and the afforestation one is 27%; as regard commission error, the afforestation registered the smallest value (14%), while the forest one is 24%.

Tab. 2 - Omission and commission error for each classification model.

Type (%)
Classification model non-forest afforestation forest
Omissionsoob non-forest / afforestation 26.2 34.4 -
afforestation / forest - 27.9 11.5
non-forest / forest 23.0 - 21.3
non-forest / afforestation / forest 36.1 39.3 23.0
Commissionsoob non-forest / afforestation 31.8 28.6 -
afforestation / forest - 13.7 23.9
non-forest / forest 21.7 - 22.6
non-forest / afforestation / forest 40.0 32.7 25.4

  Enlarge/Reduce  Open in Viewer

The large random forest OAoob values were reached thanks to the use of the hyperparameters with the largest OAcv, calculated using random search and 5-fold CV on the classification models, to identify the couple of hyperparameters which allow to reach the best results. Fig. 4 underlines that the OAcv is larger than 65% in every hyperparameters combination considered. Regarding the best iteration (value “max” in Fig. 5) the smallest OAcv is registered in the non-forest, afforestation and forest (72%) and non-forest and afforestation (75%) classification models. The model performance is larger in the classification of non-forest/forest and afforestation/forest, where the best iteration OAcv reaches 83% in the non-forest/forest classification and 87% in afforestation/forest one.

Fig. 4 - Five-fold cross validation overall accuracy (Oacv). The value “max” is referred to the best iteration. Each OAcv of every couple of hyperparameters for each classification models is displayed, grouped by accuracy.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Fig. 5 - Importance ranking of temporal predictors. The importance value (Mean Decrease Gini) of each temporal predictors in every classification model is reported.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Random forest temporal predictors importance and classification periods differences analysis

The importance ranking of temporal predictors is illustrated in Fig. 5. The most important predictors for the model afforestation/forest were the median value of the green (MDG = 2.62) and the blue (MDG = 2.62) bands, the Tasseled Cap Greenness R (MDG = 3.04), and the slope value of the Tasseled Cap Brightness (MDG = 1.76) and the green band (MDG = 1.72). The other predictors had smaller importance scores. In the non-forest and forest models, only two predictors were indicated as more important, that is the median value of the green (MDG = 10.35) and the blue (MDG = 19.36) bands; the same ranking was registered in the non-forest, afforestation, and forest models, even if the other predictors had a larger score than the ones of forest, non-forest classification. Finally, the first predictor in the importance ranking of non-forest and afforestation models is the R value of the EVI (MDG = 9.90) and, to a lesser extent, the slope value of the blue band (MDG = 2.69).

Fig. 6 displays the temporal predictors of the non-forest, afforestation and forest, related to the variables which were considered as most important in the random forest. Tasseled Cap Greeness R value (TCG - Fig. 6d), registered the largest difference between data, in particular the non-forest reaches the largest values, as in the EVI (Fig. 6c). As regard the median values of blue (Fig. 6a) and green (Fig. 6b) bands and the slope values of green (Fig. 6e) band and TCB (Fig. 6f), the forest reaches the largest values, and the afforestation and the non-forest are quite similar.

Fig. 6 - Density distribution plots. The most important temporal predictors of non-forest, afforestation and forest are displayed; plots allow to estimate the differences of the classification periods.

  Enlarge/Shrink   Download   Full Width  Open in Viewer


Afforestation monitoring is important because it influences ecosystem services but it can also cause some issues, as land fragmentation and habitat loss. In this research we aimed to demonstrate the efficacy of remote sensing data in the assessment of the difference in photosynthetic activity and in the automatic classification of non-forest, afforestation and forest areas. To classify them, we analysed 36 years of Landsat images using sixty-one training areas and a random forest algorithm. The large overall accuracy reached with the classification algorithm supports the approach chosen, confirming the efficacy of the remote sensing in the automatic classification between non-forest, afforestation and forest.

The methodology by using Landsat images allows analysing long time series of images, which is essential to assess forest evolution, which is generally slow. Moreover, both the first analysis on temporal predictors, and the random forest algorithm are easily replicable, also using different input images. To limit the overfitting due to the use of neighbouring pixels as training and validation data the data were elaborated at polygon level, considering the training data boundaries. Furthermore, random search and k-fold cross validation have been useful to identify the combination of random forest hyperparameters with the highest accuracy, to ensure the maximum performance of the algorithm. Considering the results of the random forest validation, the difference in OAoob values can be due to the similarity in photosynthetic activity in non-forest and afforestation land covers, which registered smaller accuracy values, while forest and non-forest or forest and afforestation classification models gives better results because there is a stronger difference in the spectral response of the two classes considered. It is notable that, even if the k-fold cross validation and the Out-Of-Bag accuracy results are not comparable, the OAcv and OAoob ranking is the same, with the largest accuracy value registered in distinguishing between afforestation and forest and the smallest ones in distinguishing among non-forest, afforestation and forest. This is important to support the results and the methodology, because the same result is obtained from two independent validation models.

The analysis of the temporal predictors importance allowed us to identify those variables which were more significant to separate the three classification periods. The median value of the green and the blue bands was part of the main predictors in almost every distinction considered; this probably concerned the fact that the blue band is associated with soil moisture, while the green band is associated with the reflectance value of the vegetation: the different water content and the different reflectance value are therefore two important factors to examine to classify non-forest, afforestation and forest. Considering the other predictors in the importance ranking, the Vegetation Spectral Indices and the Tasseled Cap in the first positions were those ones which are related to the blue band (as TCB, EVI), to the green band (TCG) and, to a lesser extent, those related to nir and red bands, as the Normalized Difference Vegetation Index.

These are important information to set the best combination of parameters to obtain the maximum accuracy in the distinction of non-forest, afforestation and forest, using a small number of input data. The results of the analysis of the temporal predictors confirm that the remote sensing data allow to assess the differences in the trend of photosynthetic activity in relation to non-forest, afforestation and forest. The median value of blue and green bands, r value of EVI, TCG, TCB and green slope values, which were the most important variables in the random forest ranking, showed the most significant dissimilarity in the photosynthetic activity among classes. This can be due to several reasons, not only to the spectral response of the different vegetation type, but also to their phenology; in fact, a non-forest land cover, like pastures or other grasslands, could be more influenced by the microclimate or the soil moisture than the afforestation or forest land cover. Similarly, shrublands could be more susceptible to local weather condition than the forest, which is the most stable land cover.

There are few studies on the use of remote sensing to assess afforestation areas, and most of them evaluate the forest gain by comparing different land cover maps. For example, Shen et al. ([44]) estimated afforestation areas at a regional scale using annual forest mask and combining them to extract the afforestation areas. Brovelli et al. ([7]) and Townshend et al. ([49]) assessed forest cover change using machine learning and remote sensing data through the comparison of two forest cover maps. Other researches examine the forest gain at global level using remote sensing, comparing forest and non-forest land cover, but they do not consider the afforestation land cover ([22], [28], [45]).

One of the limits of the method proposed is the medium spatial resolution of images: in fact, Landsat images resolution is 30 meters, and this could cause the omission of the smaller elements. This could be overcome by integrating other data, like radar, lidar, or other multispectral images. Another issue could occur by increasing the area of study outside the training areas; in fact, it will be necessary to consider potential commission or omission errors which could be due for example to the orography, which could reduce the accuracy because of slope and shadows.

The results of the present research offer the possibility to deepen the existing analysis approaches, directly highlighting the afforestation, a land cover class which is not stable through time, but represent a transition from non-forest to forest cover. Furthermore, it will be possible to elaborate an afforestation areas map, at national or European level. The new cartographic data could satisfy the requirements proposed by European legislation for monitoring activities, being aligned to the EIONET Action Group on Land monitoring in Europe (EAGLE) Concept. EAGLE is based on a classification system which considers the distinction between land cover and land use and it is composed of three descriptors (land cover components, land use attributes and further characteristics), combined to define a classification system suitable to be integrated to existing classes, maintaining the three components independent ([3]). The classification structure of the EAGLE matrix specifically includes the afforestation land cover class. Considering this new classification method, the map will be aligned not only to the existing cartographic data, but also to the future steps of the land cover monitoring and classification planned by Europe and Italy, being easily updated in all the descriptors.


Non-forest, afforestation and forest can be automatically predicted using Landsat data and classification models aimed to evaluate the different responses of the photosynthetic activity of non-forest, afforestation and forest areas. Since this approach is suitable for application at different scales of analysis and thanks to the results obtained in this study, the methodology will be applied to the entire Italian national area, to locate the potential afforestation areas outside the training areas. This step will allow to create a map of afforestation areas across Italy. This product follows the indications given in article 15 of the Italian national legislation on forestry and forestry supply chain (legislative decree 34/2018), which encourages the creation of a geo-referenced forest mapping and contributes advantageously to forest planning and management activities. In fact, the expansion of forest areas is often due to the agricultural and pasture surfaces loss ([34], [40]), and this could bring to an uncontrolled vegetation growth, which needs to be located and monitored, to avoid negative consequences, such as fires, that could be more destructive if the vegetation is not managed.

The afforestation area assessment can be important to analyze the forest growth trend, allowing to make forecasts on the forest evolution, both in quantitative and qualitative terms. Moreover, the integration between the map and other data will allow to study the evolution of the forest areas to assess the forest ecosystem services.


Conceptualization: SF, AC, MMa, GSM, GC, MMu; methodology: SF; statistical analysis and results elaborations: AC and SF; writing-original draft preparation: SF, AC, VF, GC, GLS, CC, LC, MMa, GC, GSM, MMu. All authors have read and agreed to the published version of the manuscript.


Agrillo E, Filipponi F, Pezzarossa A, Casella L, Smiraglia D, Orasi A, Attorre F, Taramelli A (2021). Earth observation and biodiversity big data for forest habitat types classification and mapping. Remote Sensing 13 (7): 1231.
CrossRef | Gscholar
Angelucci G (2011). Dinamica di vegetazione in aree di post-abbandono della pianura padana. Report annuale [Vegetation dynamics in post-abandonment areas of the Po Valley. Annual report]. PhD thesis, Archivio Istituzionale della Ricerca (AIR), Università di Milano, Italy, pp. 115. [in Italian]
Arnold S, Kosztra B, Banko G, Smith G, Hazeu G, Bock M (2013). The EAGLE concept - A vision of a future European Land Monitoring Framework. In: “EARSeL Symposium Proceedings 2013, Towards Horizon 2020”. Matera (Italy) 3-6 June 2013. EARSeL and CNR, Rome, Italy, pp. 551-568.
Belgiu M, Dragu L (2016). Random forest in remote sensing: a review of applications and future directions. ISPRS Journal of Photogrammetry and Remote Sensing 114: 24-31.
CrossRef | Gscholar
Bergstra J, Bengio Y (2012). Random search for hyper-parameter optimization. Journal of Machine Learning Research 13: 281-305.
Online | Gscholar
Breiman L, Cutler A (2001). Random Forest, machine learning. Statistics Department, University of California, Berkely, CA, USA, pp. 33.
Brovelli MA, Sun Y, Yordanov V (2020). Monitoring forest change in the Amazon using multi-temporal remote sensing data and machine learning classification on Google Earth Engine. ISPRS International Journal of Geo-Information 9 (10): 1-21.
CrossRef | Gscholar
Czimczik CI, Mund M, Schulze ED, Wirth C (2005). Effects of reforestation, deforestation, and afforestation on carbon storage in soils. SEB Experimental Biology Series, Taylor and Francis, Abingdon-on-Thames, UK, pp. 319-330.
CrossRef | Gscholar
De Fioravante P, Luti T, Cavalli A, Giuliani C, Dichicco P, Marchetti M, Chirici G, Congedo L, Munafò M (2021). Multispectral sentinel-2 and sar sentinel-1 integration for automatic land cover classification. Land 10 (6): 1-35.
CrossRef | Gscholar
FAO (2001). Global forest resources assessment 2000: main report. Food and Agriculture Organization, Rome, Italy.
Online | Gscholar
FAO (2020). Global forest resources assessment 2020: main report. Food and Agriculture Organization, Rome, Italy, pp. 184.
CrossRef | Gscholar
Filipponi F, Valentini E, Xuan AN, Guerra Wolf C F, Andrzejak M, Taramelli A (2018). Global MODIS fraction of green vegetation cover for monitoring abrupt and gradual vegetation changes. Remote Sensing 10 (4): 653.
CrossRef | Gscholar
Forzieri G, Girardello M, Ceccherini G, Spinoni J, Feyen L, Hartmann H, Beck PSA, Camps-Valls G, Chirici G, Mauri A, Cescatti A (2021). Emergent vulnerability to climate-driven disturbances in European forests. Nature Communications 12 (1): 1-12.
CrossRef | Gscholar
Francini S (2021). BAP-GEE - A Google Earth Engine application for Best Available Pixel composites calculation, visualization, calibration, and download. Web site.
Online | Gscholar
Francini S, McRoberts RE, Giannetti F, Mencucci M, Marchetti M, Scarascia Mugnozza G, Chirici G (2020). Near-real time forest change detection using PlanetScope imagery. European Journal of Remote Sensing 53 (1): 233-244.
CrossRef | Gscholar
Francini S, Hermosilla T, Coops N, White J, Wulder M, Chirici G (2021a). A Google Earth Engine application for Best Available Pixel composites calculation, visualization, calibration, and download. Web site.
Online | Gscholar
Francini S, McRoberts RE, Giannetti F, Marchetti M, Scarascia Mugnozza G, Chirici G (2021b). The Three Indices Three Dimensions (3I3D) algorithm: a new method for forest disturbance mapping and area estimation based on optical remotely sensed imagery. International Journal of Remote Sensing 42 (12): 4697-4715.
CrossRef | Gscholar
Giannetti F, Pegna R, Francini S, McRoberts RE, Travaglini D, Marchetti M, Mugnozza GS, Chirici G (2020). A new method for automated clearcut disturbance detection in Mediterranean coppice forests using Landsat time series. Remote Sensing 12 (22): 1-23.
CrossRef | Gscholar
Gorelick N, Hancher M, Dixon M, Ilyuschenko S, Thau D, Moore R (2017). Google Earth Engine: planetary-scale geospatial analysis for everyone. Remote Sensing of Environment 2020: 18-27.
CrossRef | Gscholar
Griffiths P, Van Der Linden S, Kuemmerle T, Hostert P (2013). Erratum: a pixel-based landsat compositing algorithm for large area land cover mapping. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing 6 (5): 2088-2101.
CrossRef | Gscholar
Gálos B, Hagemann S, Hänsler A, Kindermann G, Rechid D, Sieck K, Teichmann C, Jacob D (2013). Case study for the assessment of the biogeophysical effects of a potential afforestation in Europe. Carbon Balance and Management 8 (1): 1-12.
CrossRef | Gscholar
Hansen MC, Potapov P, Moore R, Hancher M, Turubanova SA, Tyukavina A, Thau D, Stehman S, Goetz SJ, Loveland TR, Kommareddy A, Egorov A, Chini L, Justice CO, Townshend JRG (2013). High-resolution global maps of 21st-century forest cover change. Science 342 (6160): 850-853.
CrossRef | Gscholar
Hawrylo P, Francini S, Chirici G, Giannetti F, Parkitna K, Krok G, Mitelsztedt K, Lisanczuk M, Sterenczak K, Ciesielski M, Wezyk P, Socha J (2020). The use of remotely sensed data and Polish NFI plots for prediction of growing stock volume using different predictive methods. Remote Sensing 12 (20): 1-20.
CrossRef | Gscholar
Hermosilla T, Wulder MA, White JC, Coops NC, Hobart GW (2015a). An integrated Landsat time series protocol for change detection and generation of annual gap-free surface reflectance composites. Remote Sensing of Environment 158: 220-234.
CrossRef | Gscholar
Hermosilla T, Wulder MA, White JC, Coops NC, Hobart GW (2015b). Regional detection, characterization, and attribution of annual forest change from 1984 to 2012 using Landsat-derived time-series metrics. Remote Sensing of Environment 170: 121-132.
CrossRef | Gscholar
Huang H, Chen Y, Clinton N, Wang J, Wang X, Liu C, Gong P, Yang J, Bai Y, Zheng Y, Zhu Z (2017). Mapping major land cover dynamics in Beijing using all Landsat images in Google Earth Engine. Remote Sensing of Environment 202: 166-176.
CrossRef | Gscholar
Jönsson P, Cai Z, Melaas E, Friedl MA, Eklundh L (2018). A method for robust estimation of vegetation seasonality from Landsat and Sentinel-2 time series data. Remote Sensing 10 (4): 635.
CrossRef | Gscholar
Kim DH, Sexton JO, Noojipady P, Huang C, Anand A, Channan S, Feng M, Townshend JR (2014). Global, Landsat-based forest-cover change from 1990 to 2000. Remote Sensing of Environment 155: 178-193.
CrossRef | Gscholar
Koch B (1999). The contribution of remote sensing for afforestation and forest biodiversity assessment. In: Proceedings of the conference “Remote Sensing and Forest Monitoring”. Rogow (Poland) 1-3 June 1999. T. Zawila-Niedzwiecki, M. Brach, Poland, pp. 3-13.
Laurin GV, Francini S, Luti T, Chirici G, Pirotti F, Papale D (2021). Satellite open data to monitor forest damage caused by extreme climate-induced events: a case study of the Vaia storm in Northern Italy. Forestry 94 (3): 407-416.
CrossRef | Gscholar
Luti T, De Fioravante P, Marinosci I, Strollo A, Riitano N, Falanga V, Mariani L, Congedo L, Munafò M (2021). Land consumption monitoring with SAR data and multispectral indices. Remote Sensing 13 (8): 1-26.
CrossRef | Gscholar
Marchetti M, Bertani R, Corona P, Valentini R (2012). Changes of forest coverage and land uses as assessed by the inventory of land uses in Italy. Forest@ - Journal of Silviculture and Forest Ecology 9 (4): 170-184. [in Italian with English summary]
CrossRef | Gscholar
Mather AS (2000). Afforestation: progress, trends and policies. In: Proceedings of the “Scientific Symposium”. Freiburg (Germany) 26-17 February 2000. European Forest Institute Proceedings no. 35, pp. 11-19.
Munafò M (2018). Territorio: processi e trasformazioni in Italia [Territory: processes and transformations in Italy]. Report 296/2018, Istituto Superiore per la Protezione e la ricerca Ambientale (ISPRA), Sistema Nazionale per la Protezione dell’Ambiente (SNPA), Rome, Itay, pp. 88. [in Italian]
Nguyen HTT, Doan TM, Tomppo E, McRoberts RE (2020). Land use / land cover mapping using multitemporal Sentinel-2 imagery and four classification. Remote Sensing 12: 1-27.
Nicodemus KK (2011). Letter to the editor: on the stability and ranking of predictors from random forest variable importance measures. Briefings in Bioinformatics 12 (4): 369-373.
CrossRef | Gscholar
Palmero-Iniesta M, Espelta JM, Gordillo J, Pino J (2020). Changes in forest landscape patterns resulting from recent afforestation in Europe (1990-2012): defragmentation of pre-existing forest versus new patch proliferation. Annals of Forest Science 77 (2): 971.
CrossRef | Gscholar
Persson M, Lindberg E, Reese H (2018). Tree species classification with multi-temporal Sentinel-2 data. Remote Sensing 10 (11): 1-17.
CrossRef | Gscholar
Qiu B, Zou F, Chen C, Tang Z, Zhong J, Yan X (2018). Automatic mapping afforestation, cropland reclamation and variations in cropping intensity in central east China during 2001-2016. Ecological Indicators 91: 490-502.
CrossRef | Gscholar
RAF (2019). RaFITALIA 2017-2018. Rapporto sullo stato delle foreste e del settore forestale in Italia [Report on the state of forests and the forestry sector in Italy]. Compagnia delle Foreste, Arezzo, Italy, pp. 284. [in Italian]
Raha AK, Mishra AV, Das S, Zaman S, Ghatak S, Bhattacharjeef S, Raha S, Mitra A (2010). Time series analysis of forest and tree cover of West Bengal from 1988 to 2010, using RS/GIS, for monitoring afforestation programmes. The Journal of Ecology 108: 255-265.
Sallustio L, Munafò M, Riitano N, Lasserre B, Fattorini L, Marchetti M (2016). Integration of land use and land cover inventories for landscape management and planning in Italy. Environmental Monitoring and Assessment 188 (1): 1-20.
CrossRef | Gscholar
Schmidt G, Jenkerson C, Masek J, Vermote E, Gao F (2013). Landsat ecosystem disturbance adaptive processing system (LEDAPS) algorithm description. Report 2013-1057, USGS Science for a Changing World, US Geological Survey, Washington, DC, USA, pp. 17.
Shen W, Li M, Huang C, Tao X, Li S, Wei A (2019). Mapping annual forest change due to afforestation in Guangdong Province of China using active and passive remote sensing data. Remote Sensing 11 (5): 1-21.
CrossRef | Gscholar
Shimada M, Itoh T, Motooka T, Watanabe M, Shiraishi T, Thapa R, Lucas R (2014). New global forest/non-forest maps from ALOS PALSAR data (2007-2010). Remote Sensing of Environment 155: 13-31.
CrossRef | Gscholar
Smiraglia D, Filipponi F, Mandrone S, Tornato A, Taramelli A (2020). Agreement index for burned area mapping: integration of multiple spectral indices using Sentinel-2 satellite images. Remote Sensing 12 (11): 1862.
CrossRef | Gscholar
Spadoni GL, Cavalli A, Congedo L, Munafò M (2020). Analysis of Normalized Difference Vegetation Index (NDVI) multi-temporal series for the production of forest cartography. Remote Sensing Applications: Society and Environment 20: 100419.
CrossRef | Gscholar
Spinsi A, Botarelli L, Marletto V (2012). Indici vegetazionali da satellite per il monitoraggio in continuo del territorio [Satellite-based vegetation indices for continuous land monitoring]. Italian Journal of Agrometeorology 3: 49-55. [in Italian]
Townshend JR, Masek JG, Huang C, Vermote EF, Gao F, Channan S, Sexton JO, Feng M, Narasimhan R, Kim D, Song K, Song D, Song XP, Noojipady P, Tan B, Hansen MC, Li M, Wolfe RE (2012). Global characterization and monitoring of forest cover using Landsat data: opportunities and challenges. International Journal of Digital Earth 5 (5): 373-397.
CrossRef | Gscholar
UN-FCCC (2013). Afforestation and reforestation projects under the clean development mechanism. A reference manual. United Nations Framework Convention on Climate Change, Bonn, Germany, pp. 3-63.
White JC, Wulder MA, Hermosilla T, Coops NC, Hobart GW (2017). A nationwide annual characterization of 25 years of forest disturbance and recovery for Canada using Landsat time series. Remote Sensing of Environment 194: 303-321.
CrossRef | Gscholar
Whittaker G, Confesor R, Di Luzio M, Arnold JG (2010). Detection of overparameterization and overfitting in an automatic calibration of SWAT. Transactions of the ASABE 53 (5): 1487-1499.
CrossRef | Gscholar

Authors’ Affiliation

Alice Cavalli 0000-0002-5460-1245
Mauro Maesano 0000-0002-4325-951X
Giuseppe Scarascia Mugnozza 0000-0003-0357-4360
Department of Innovation in Biology, Agri-Food and Forest systems - DIBAF, University of Tuscia, v. San Camillo de’ Lellis snc, I-01100 Viterbo (Italy)
Saverio Francini 0000-0001-6991-0289
Gherardo Chirici 0000-0002-0669-5726
Fondazione per il futuro delle città, Firenze (Italy)
Saverio Francini 0000-0001-6991-0289
Claudia Cocozza 0000-0002-0167-8863
Gherardo Chirici 0000-0002-0669-5726
Department of Agricultural, Food and Forestry Systems, University of Florence (Italy)
Giulia Cecili 0000-0002-8199-7660
Valentina Falanga 0000-0003-0454-8850
Gian Luca Spadoni 0000-0001-6083-6051
Dept. of Biosciences and Territory, University of Molise, c.da Fonte Lappone, I-86090 Pesche, IS (Italy)
Luca Congedo 0000-0001-7283-116X
Michele Munafò 0000-0002-3415-6105
Italian Institute for Environmental Protection and Research - ISPRA, v. Vitaliano Brancati 48, I-00144 Rome (Italy)

Corresponding author

Saverio Francini


Cavalli A, Francini S, Cecili G, Cocozza C, Congedo L, Falanga V, Spadoni GL, Maesano M, Munafò M, Chirici G, Scarascia Mugnozza G (2022). Afforestation monitoring through automatic analysis of 36-years Landsat Best Available Composites. iForest 15: 220-228. - doi: 10.3832/ifor4043-015

Academic Editor

Agostino Ferrara

Paper history

Received: Dec 20, 2021
Accepted: May 06, 2022

First online: Jul 12, 2022
Publication Date: Aug 31, 2022
Publication Time: 2.23 months

© SISEF - The Italian Society of Silviculture and Forest Ecology 2022

  Open Access

This article is distributed under the terms of the Creative Commons Attribution-Non Commercial 4.0 International (https://creativecommons.org/licenses/by-nc/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Creative Commons Licence

Breakdown by View Type

(Waiting for server response...)

Article Usage

Total Article Views: 234
(from publication date up to now)

Breakdown by View Type
HTML Page Views: 0
Abstract Page Views: 0
PDF Downloads: 224
Citation/Reference Downloads: 3
XML Downloads: 7

Web Metrics
Days since publication: 149
Overall contacts: 234
Avg. contacts per week: 10.99

Article citations are based on data periodically collected from the Clarivate Web of Science web site
(last update: Nov 2020)

(No citations were found up to date. Please come back later)


Publication Metrics

by Dimensions ©

List of the papers citing this article based on CrossRef Cited-by.


iForest Similar Articles


This website uses cookies to ensure you get the best experience on our website. More info