Harmonized Classiﬁcation of Forest Types in the Iberian Peninsula Based on National Forest Inventories

: National Forest Inventories (NFIs) collect and provide a large amount of information regarding the forest volume, carbon stocks, vitality, biodiversity, non-wood forest products and their changes. Forest stands variables data are paramount to understanding their composition, especially on those related with understory characteristics and the coverage of species according to canopy layers; they are essential to assess biodiversity and to support forest management. At the same time, these inventories allow the development of harmonized forest descriptions beyond the national scale. This study aims to develop a homogeneous characterization of the Iberian Peninsula’s forests, in order to classify and identify the forest types. For this purpose, harmonized data from NFIs of Portugal and Spain were used to assess the composition of species, dominance and the percentage of cover for each species in a vertical space deﬁned by seven canopy layers. Using the “K-means” clustering algorithm, a set of clusters was identiﬁed and georeferenced using forest polygons from land use and cover maps of both countries. The interpretation and description of the clusters lead to the establishment of 28 forest types that characterize all of the Iberian Peninsula forests. Each forest area has been described through one of the forest types and their relation with other ecological characteristics of the stands was analyzed. Shrubs formations are generally widely distributed in the forest area of the Iberian Peninsula, however their abundance in terms of cover is lower in comparison with tree species. Around 71% of the forest types are dominated by trees, mainly species from the genera Pinus and Quercus , and 21% are dominated by shrub formations with species of Ulex spp., Cytisus spp., and Cistus spp. The Quercus ilex s.l. L. and Pinus pinaster Aiton are the common species of importance for both NFIs. The results represent a powerful and homogenous multi-use tool describing the Iberian Peninsula’s forestlands with applications on landscape analysis, forest management and conservation. This information can be used for comparisons at larger scales, allowing cross-border analysis in relation to various aspects, such as hazards and wildﬁres, as well as management and conservation of forest biodiversity. The developed method is adaptable to an updated dataset from more recent NFIs and to other study areas.


Introduction
National Forest Inventories (NFIs) are implemented in order to provide information about the status of forest resources at the sub-national and national levels and to assist not only in the formulation of forestry policies, but also in the strategy for industrial policies. These instruments constitute the main source of forest information in terms of exhaustiveness and precision in each country, providing information on the characteristics of the forest stands in a larger area. Over time, the NFIs have evolved gradually and changed their assessment from forest volume information to sustainable ecosystem indicators, which include new variables to estimate carbon and biodiversity that support forest management and the implementation of fuel management programs [1]. Recent NFIs include measurements concerning species composition on different canopy layers and understory characteristics, which are crucial to characterize the forest structure and its functions on the ecosystems [2].
Inventory cycles, with sample-based techniques, are currently used across European and American countries. However, due to the different historical background, countries use their own definitions, methodologies and variables thresholds to fulfil their national forest policies needs and forest programs. In Europe, with the establishment of the European National Forest Inventory Network (ENFIN) in 2003, and later the COST Action E43, work has been done to harmonize estimations, develop methods, concepts and definitions to be employed in the NFIs from different countries in Europe so the information would become comparable [1]. Harmonization is the process to achieve reference definitions and methods to produce comparable estimates and acknowledges that individual countries have developed the unique features of their NFIs for specific purposes and are justified in their desire to maintain them [1,3,4]. The harmonization process implies the use of a bridge function or conversion from local to the established reference definition [4,5]. The harmonized information support international reporting requirements, which have become more frequent due to the increasing agreements and commitments such as the FAO Forest Resources Assessment and the UNECE Temperate and Boreal Forest Recourses Assessment, and support also political decision making. At the European scale, several harmonization studies have been developed on forest area, growing stock, biodiversity indicators and conservation status taking into account data from NFIs [2,[6][7][8][9].
The information related to forest fires prevention at large scale has been subject to a wider interest lately and, likewise, the European Forest Fire Information System (EFFIS) has been able to gather together most of the European countries to produce harmonized reporting on forest fire information [6]. This harmonized information is fundamental to having cross-border and integrated data about wildfires that can potentially assist in firefighting.
Even if Portugal and Spain are quite different in size and population, and the national forest monitoring processes for each country have been influenced by different historical, cultural, political and economic contexts, forests in the Iberian Peninsula share similar environmental conditions, fire propensity, stand structure characteristics and common species [10].
Wildfires associated to rural abandonment represent one of the major concerns for the future of the forest and forestry sector in the Iberian Peninsula. The burned area of forest stands and shrubs was about 442 thousand hectares in Portugal [11] and 178 thousand hectares in Spain [12], in the year 2017, which represents around 14% and 1% of the forest areas, respectively. Thus, the analysis of the forest type classification and the wildfire probability for the Iberian Peninsula needs to be treated as a whole.
Iberian Peninsula forests are divided into the predominant Mediterranean and Atlantic bioregions [13]. Besides these two regions there is a smaller alpine region forest area. The Atlantic part covers a continuous area in the northern and north-western territories of the Peninsula, with larger areas of heathlands, oak forests and grasslands in the mountains with scattered juniper on rocky areas among the pine forests. The Mediterranean Iberia is dominated by sclerophyllous, semi-deciduous forests and shrublands, and the mountainous areas occupied by different species of pine trees. As a result of human impact, evergreen broadleaf forest have been converted to shrubland areas or to silvopastoral oak woodlands in the western part of the Iberian Peninsula [14,15]. In Portugal, forests cover 3.2 million hectares (35.3% of total land) and other wooded land 1.7 million hectares [16], whereas in Spain forests cover 18.4 million (36.9% of total land) and other wooded land 9.2 million hectares. In Spain there is a target to increase the forest area by 4 million hectares by 2032 [16]. The FAO definition of forest is used for both countries [17]. The forest area available for wood supply in 2015 was 2.1 million hectares in Portugal and 14.7 million hectares in Spain.
The information on NFIs can also contribute to characterize the fuel present in each vegetation type based on the understory fuel structure. Also, the forest characteristics related to the vertical and horizontal structure of the stands are one key to predict the behaviour of a wildfire. The use of NFIs data in relation to fires has been explored by several authors [7,[18][19][20][21]. The horizontal structure is often assessed by aerial photographs or satellite information and NFI plots can be used to complement and validate those assessments [22,23]. However, in spite of recent developments with LIDAR, information on the vertical structure of forests is obtained or validated mainly from NFIs. Apart from crown characteristics, the shrub cover also constitutes an important ecological indicator, as suggested by multiple international reporting organizations (e.g., United Nations Framework Convention on Climate Change), but the integration of the information from the upper strata in the same analysis for a large-scale reporting is generally lacking. Non-tree vegetation species, e.g., shrubs, play a great ecological role on biodiversity as bioindicators for erosion and landscape fire recovery, and therefore are important factors in forest management [24,25]. It is essential to assess the whole vertical structure of the forests [2,26] and the information provided by the NFIs is very useful to represent these structures on larger areas.
The characteristics and location of the forest stands are essential to classify ecological forest types and the information in the NFIs is the only currently available regarding forest resources and distribution at country and sub-country levels. Forest classification defined by its composition and/or site factors, support silvicultural and forest ecosystem planning decisions from local to global scales [27]. Despite the NFIs being an extremely rich source of information about forests with different potential applications [26], data about composition of the canopy layers is collected in a non-harmonized way in the European countries, which compromises the comparability and reinforces the need to develop harmonized methodologies on this matter to assess ecological forest types and information on fuel loads. Harmonized analysis of a continuous landscape is essential for comparison matters in a larger scale [28]. Portugal and Spain have available information, collected on their NFIs, on cover, species composition, both understory and canopy levels, and data on these matters are relevant to describe shrub communities and classify the vegetation in forest types. However, their information, definitions and rules are country-specific, each country has their own method to monitor species composition and structure [29].
To address this situation, we implemented a harmonization process using information from the Portuguese and Spanish NFIs. This methodology can be easily applied to other European countries. The aim of the study was threefold: (i) to present an approach based on NFIs data of Portugal and Spain to harmonize the data of stands composition coverage according to their height classes; (ii) to identify the forest types in the Iberian Peninsula using harmonized data; (iii) to analyze the relation of the forest types with other ecological characteristics of the stands like geographical characteristics and wildfire probability.

Steps of the Harmonization Process
The harmonized methodological process was carried out within the framework of the DIABOLO project. We took the methodology developed within the COST Action 43 that comprises four steps [30] and applied it in the Iberian Peninsula. In the first step, we studied the different international and national definitions related with forest and shrub types, aiming to define and determine common forest and shrubs types. To achieve this goal, we assessed the information available both in the Portuguese and Spanish National Forest Inventories (NFIs) in terms of type of species, tree and stand level information, composition of the species through the different canopy layers (vertical structure), thresholds of each canopy and the distribution of main forest across countries. An analytical decomposition of each definition was then made to be harmonized between countries. In the second step, we established the reference definition for the information mentioned on step one. In the third step, we built bridges between the national definitions of vertical structure in each NFI. The reference definition to transform national data into comparable data was made in the fourth step. In the following sections, we present the description of data sources and the methods applied on the harmonization steps.

NFI Data for Portugal and Spain
NFIs in the Iberian Peninsula have been producing statistics and cartography based on the abundance, state and condition of national forests resources [31,32]. The Portuguese National Forest Inventory (PTNFI) was started in 1965 and covers the entire territory of mainland Portugal. The information has been updated with a periodicity of approximately 10 years. So far, six NFIs were completed. The Spanish National Forest Inventory (SNFI) covers all forest land in Spain, was started also in the 1965 and the last NFI, the fourth, is currently on-going.
Data for this study were obtained from the Fifth Portuguese National Forest Inventory (PTNFI5, 2005(PTNFI5, -2006 and the Third Spanish National Forest Inventory (SNFI3, 1997(SNFI3, -2007. An overview about the NFI sampling methods relevant to characterize those NFIs is given in Table 1. In the PTNFI5 shrub attributes are measured in 10 m circular radius sample plots. The SNFI3 visually assesses mean height and cover for each shrub species using a percentage scale with 1% interval width [2]. In the SNFI3 data collection can be categorized in three main classes: stand description, measured tree data and biodiversity assessment. Saplings and regeneration are recorded in the 5 m radius subplots. All the continental areas of Portugal and Spain with more than 0.5 hectares with trees higher than 5 m and a canopy cover of more than 10 percent, consistent with FAO definition of Forest Area [17], were identified and analyzed and constitute the study area. The database used for the harmonization process has a total of 7935 inventory plots from Portugal and 60,607 inventory plots from Spain. Forest areas in the Macaronesia Region (Canary, Azores and Madeira islands) and the Mediterranean Balearic islands were not considered.
A common legend of species or group of species in the Iberian Peninsula has been obtained, matching and aggregating, through taxonomic and ecological criteria, tree and shrub species/formations  (Table S1). The nomenclature of species is in accordance with Flora Iberica [33].
Shrub and grass species of PTNFI [34] are consistent with the SNFI3 and Spanish Forest Map (SFM25) [35] hierarchical categorization. The matching of species and formation from both NFIs led to the definition of 14 grass and shrub formation identified by codes according to SNFI3 and SFM25 species categorization (Table S2).

Composition of Stands
In both inventories, PTNFI5 and SNFI3, the following data were recorded for each sample tree: species, diameter at breast height, height, abundance as well information regarding the composition of the species through the different canopy layers, the so-called vertical structure. Countries have different thresholds of the vertical structure height classes, so we propose a harmonization process to establish, by height classes, the percentage of cover by species, tree and shrubs. We followed the description adopted in PTNFI, consisting in seven height strata [36] from the highest (stratum 1: >16 m) to the lowest (stratum 7: 0-0.5 m) ( Figure 1). These thresholds were already used in Portugal and were adapted for the Spanish NFI data ( osition of Stands th inventories, PTNFI5 and SNFI3, the following data were recorded for each sam iameter at breast height, height, abundance as well information regarding the com ecies through the different canopy layers, the so-called vertical structure. Count thresholds of the vertical structure height classes, so we propose a harmonizatio sh, by height classes, the percentage of cover by species, tree and shrubs. We foll n adopted in PTNFI, consisting in seven height strata [36] from the highest (stratu lowest (stratum 7: 0-0.5 m) ( Figure 1). These thresholds were already used in Port pted for the Spanish NFI data (  In the PTNFI5, vertical structure was assessed by the cover percentages for each of the three most abundant species (tree or shrub) per height class. In the SNFI3, cover percentage per layer was not estimated directly. We developed models to estimate the geometry of the tree crowns, crown diameter (CD, m) and crown length (CL, m), of each tree using information from four sample trees systematically selected on each sample plot of the SNFI2 to complete a database with more than 255,000 sample trees covering the distribution area of all the forest species and all the combinations of age, stand density and site qualities. The estimates of these variables were used to classify the trees in each one of the seven height strata. Moreover, mean height and cover for each shrub species using a percentage scale with 1% interval width were recorded in the 10 m radius subplots in the SNFI3 [2]. These measures were used to classify the shrubs species in each stratum.
For the adjustment of the models, we included a variable related to the density of the stand, as competition for space limits the growth of individual canopies and it is well known that crown diameter and crown length decreases when stands are dense. In fact, if the distribution is completely regular in a square grid, each tree would have an available square area to grow that is the inverse of the tree density, and the square root of that area is the distance between each tree and its four closest neighbours ( Figure 2). Equation (1) expresses that relationship: where DIST = reference distance (m), NHA = number of trees (trees·ha −1 ).
Forests 2020, 11, x FOR PEER REVIEW 6 of 22 In the PTNFI5, vertical structure was assessed by the cover percentages for each of the three most abundant species (tree or shrub) per height class. In the SNFI3, cover percentage per layer was not estimated directly. We developed models to estimate the geometry of the tree crowns, crown diameter (CD, m) and crown length (CL, m), of each tree using information from four sample trees systematically selected on each sample plot of the SNFI2 to complete a database with more than 255,000 sample trees covering the distribution area of all the forest species and all the combinations of age, stand density and site qualities. The estimates of these variables were used to classify the trees in each one of the seven height strata. Moreover, mean height and cover for each shrub species using a percentage scale with 1% interval width were recorded in the 10 m radius subplots in the SNFI3 [2]. These measures were used to classify the shrubs species in each stratum.
For the adjustment of the models, we included a variable related to the density of the stand, as competition for space limits the growth of individual canopies and it is well known that crown diameter and crown length decreases when stands are dense. In fact, if the distribution is completely regular in a square grid, each tree would have an available square area to grow that is the inverse of the tree density, and the square root of that area is the distance between each tree and its four closest neighbours ( Figure 2). Equation (1) expresses that relationship: where DIST = reference distance (m), NHA = number of trees (trees.ha −1 ).

Figure 2.
Stand geometry based on the mean distance between trees considering a squared grid.
For the estimation of crown length, a logistic model was fitted for each one of the species or group of species, considering the tree height as an asymptotic value and using the reference distance as independent variable (Equation (2)). It is clear from the Equation (2) that when tree density increases the value of the reference distance tends to zero and the crown ratio, i.e., the fraction of crown length and total height, tends to the expression on Equation (3): where CL = crown length (m), HT = tree height (m), DIST = reference distance (m), = crown ration, a0 and a1 = parameters to be estimated.
Parameter a0 can be interpreted as a coefficient related to the minimum crown ratio within a canopy. Small values of a0 indicate deep canopies. Parameter a1 represents the effect of competition and it is associated with a negative sign as when the distance between trees increases the denominator of the equation decreases resulting in the increase of the estimated crown length for a given tree height. Smaller values for a1 (closer to zero) indicate smaller effects of distance to neighbour trees. A value of zero for a1 would indicate that competition with neighbour trees has no effect on crown length. For the estimation of crown length, a logistic model was fitted for each one of the species or group of species, considering the tree height as an asymptotic value and using the reference distance as independent variable (Equation (2)). It is clear from the Equation (2) that when tree density increases the value of the reference distance tends to zero and the crown ratio, i.e., the fraction of crown length and total height, tends to the expression on Equation (3): where CL = crown length (m), HT = tree height (m), DIST = reference distance (m), CL HT = crown ration, a 0 and a 1 = parameters to be estimated.
Parameter a 0 can be interpreted as a coefficient related to the minimum crown ratio within a canopy. Small values of a 0 indicate deep canopies. Parameter a 1 represents the effect of competition and it is associated with a negative sign as when the distance between trees increases the denominator of the equation decreases resulting in the increase of the estimated crown length for a given tree height. Smaller values for a 1 (closer to zero) indicate smaller effects of distance to neighbour trees. A value of zero for a 1 would indicate that competition with neighbour trees has no effect on crown length.
For the estimation of crown diameter, we used a similar approach. In this case, as we were dealing with a variable related to the horizontal structure, the tree variable used to estimate crown diameter was the diameter at breast height (DBH) of the stem, as this is the most common measure taken in forest inventories [1]. The effect of competition by neighbour trees was included based on the reference distance. The expression of the equation finally fitted was as follows: where CD = crown diameter (m), DBH = diameter at breast height (cm), DIST = reference distance (m), b 0 , b 1 and b 2 = parameters to be estimated. Parameter b 0 represents the tendency for the crown of the species to extend horizontally. In mathematical terms, it would be the estimated crown diameter (m) of a tree with 1 cm of DBH and 1 m to the nearest neighbouring trees if spacing was in a completely regular square grid ( Figure 2). Parameters b 1 and b 2 are both positive as they are related with the positive relationship between stem and crown diameters and between spacing between trees and crown diameter. Tables S4 and S5 present the characteristics of the sample trees of each species used to fit the models 2 and 4, respectively. The parameter estimates and the coefficient of determination of fitting Equations (2) and (4) to each species or groups of species are presented in Tables S6 and S7.
To describe the composition and to better understand the effect on the diversity of vertical structure, the percentages cover of forest species according to the seven height classes and litter were used as variables. Each species according to height class was considered as pseudo-species, allowing for an estimate of the proportion coverage that each species contributes to a specific height class.
The harmonized dataset has a total of 55,003 plots with information for vertical structure on the seven vertical layers (Figure 1), the shrub and/or formations and the trees and groups of trees (Tables S1 and S2), which define 217 layer/species (pseudo-species), plus one variable concerning the percentage of litter and two geographical coordinates (WGS 84), which makes a total of 220 variables. We have selected just the plots having at least a 10% of coverage value in one of the upper layers (h1, h2 and h3; h > 8 m). Even though the tree height threshold in FAO definition of forest is 5 m, an under estimation given by using an 8 m threshold has been preferred to the over estimation eventually given by a threshold at 4 m.

Forest Types and Cartography
Forest types were defined by a K-means clustering algorithm and then georeferenced using forest polygons in both countries. The forest polygons/areas were defined using land use and cover maps for Portugal "Carta de Uso e Ocupação de Solo de Portugal Continental COS" [37] and for Spain "Mapa Forestal de España" [35], both from 2007 and available as shapefiles [38]. The forest polygons were identified in both land use maps according to the FAO definition of forest area. Figure 3 presents the workflow of the different steps to establish the forest types. Specifically, we used as features for the K-means clustering the coverage of each species for each of the seven height layers per plot (pseudo-species) and also the percentage of litter. The litter information was included at the layer located on the ground level (layer seven h7). The K-means algorithm was implemented through the software SPSS 23 (IBM Analytics 2015). The value of K (number of clusters) was chosen in order to have a good distribution of sampling plots for each cluster, avoiding too many small (containing <1% of plots) or big (>10% of plots) clusters.  Each cluster is described by its centroid which is defined in a space of N dimensions corresponding to the number of variables (coverage of species by height classes and litter percentage) used for the analysis. Clusters centroids describe the average values for each variable for each forest type. A first analysis allowed us to give an interpretable name to each cluster for better identification and GIS visualizations using information about canopy cover, composition and dominant species/formation: • canopy cover: a total cover value was obtained for each type as the sum of coverage percentage for all layer/species variables. The range of values were divided in three equal groups in order to attribute the following density classes: closed forests (c), average density (no code), open forests (o); • composition: the proportion of shrubs and tree species was calculated in each forest type. The obtained range of values were divided in three equal groups which were used to name the model in relation to their dominance: shrubs (S), mixed (M) or tree dominance (T); • dominant species/formation: a code according to Tables S1 and S2 was defined taking into account the higher coverage of the species, or group of species in the cluster. In those cases where two species had a similar total value of coverage, both species were indicated. In other cases, where there is no defined dominance of one or two species, the code mix was applied. Each cluster is described by its centroid which is defined in a space of N dimensions corresponding to the number of variables (coverage of species by height classes and litter percentage) used for the analysis. Clusters centroids describe the average values for each variable for each forest type. A first analysis allowed us to give an interpretable name to each cluster for better identification and GIS visualizations using information about canopy cover, composition and dominant species/formation: dominant species/formation: a code according to Tables S1 and S2 was defined taking into account the higher coverage of the species, or group of species in the cluster. In those cases where two species had a similar total value of coverage, both species were indicated. In other cases, where there is no defined dominance of one or two species, the code mix was applied.
With the centroids analysis, it was possible to elaborate the graphs describing the composition of each forest type. The percentage of cover of each species or group of species according to each canopy layer was considered as a pseudo-species. To obtain a clearer visualization, we associated as "others" all the species which do not reach the 3% of coverage in at least one of the seven vertical height strata. The dataset was updated with the code information of each forest type, which was used to classify the plots with different colours. Using geographical coordinates included in the dataset and a generic background map, it was possible to visualize the position of each plot in the study area. The categorization of the forest types is not exhaustive and does not allow us to measure the area of extension of each forest type. A good solution to obtain a more correct and complete georeferentiation was found observing land use maps of Portugal and Spain [35,37]. These maps have the peculiarity to distinguish forest area by vegetation characteristics, including structure and species composition. This means that each polygon corresponds only to one forest type. Using a "spatial join" algorithm implemented in ArcGIS (ESRI), each forest polygon was characterized by the forest types defined for nearest NFI plot. Adjacent polygons corresponding to the same forest type, have been "dissolved" to delete unnecessary borders in the final map.

Ecological Preferences
We have analyzed the relation between each forest type with characteristics of the plots, namely latitude, elevation, slope position and slope aspect. We used geographical information from various sources to derive the values of each variable per plot. The analysis was made for dataset of the plot information. Geographical information was retrieved from the European Digital Elevation Models (EU-DEM v1.1) provided by Copernicus. Four raster files were obtained (E20N20, E30N20, E20N10, E30N10) and compiled into one to have a single DEM. Slope aspect and slope position, in raster format, was calculated using GDAL tools implemented in software QGIS. In order to have for each inventory plot the needed data, we used the plugin of QGIS named Point Sampling Tool, and added the values from the three rasters (DEM, slope aspect and slope position) to the points layer containing information of forest types (see the workflow at Figure 3). Latitude data was provided by the y-coordinate of each plot. Each geographical characteristic was arranged into the following classes: For each variable, the average for all plots of each forest type were compared with the total average value for all plots.
We also analyzed the relation between the forest types and wildfire probability for Portugal and Spain. This was established taking into account the fire dataset retrieved from the European Fire Database provided by the EFFIS. The forest fire information for the period 2005 to 2015 was derived from the daily processing of MODIS satellite imagery at 250 m ground spatial resolution, where burnt scars of approximately 30 ha in size were mapped. The identification of the plots that were later burned was done by intersecting the burned areas with the NFI plots of Portugal and Spain, using QGIS. A binary categorical variable was then created, equalling 0 if the plot was not burned later and 1 if the plot was burned later at least one time. Based on the wildfire characteristics in the Iberian Peninsula and also in the knowledge from field experts on wildfire behaviour, we assumed that all the plots that were located inside a burned area were burned [7].

Forest Types in the Iberian Peninsula
The result of the K-means analysis was the average percentage cover of species according to height class for each cluster, thus outlining the profiles of resulting forest types. The K-means clustering was made initially with a value of K = 21. Although, among the 21 clusters, one contains less than 1% of the plots and one cluster contains more than 10% of plots (12,749 plots, which represents 23% of the plots). We decided to sub-divide the bigger cluster by K = 8, obtaining a final partitioning of 28 clusters in total. These clusters characterize the Peninsula Iberian forests, with each cluster corresponding to a particular forest type. These forest types were georeferenced using forest polygons detected in land use and cover maps of Portugal and Spain [35,37]. As a result of the analysis of cluster centres, forest types have been named according to the following criteria: canopy cover, composition and dominant species or shrub formation. The new eight clusters resulting from the subdivision allowed a more specific characterization of the following vegetation presented in each forest type: Open Quercus ilex forests with low understory Each cluster contains final polygons and punctual visualization of vegetation structure. The area corresponding to each forest type was calculated, as well the percentage of the forest area represented from the total of the Iberian Peninsula and the percentage of the NFIs plots within each cluster ( Table 2). All the forest areas are in accordance with the FAO definition for forest area [17]. Because we decided not to join attributes between points and polygons at a distance greater than 10 km, this resulted in an area of about 700 km 2 , around 0.3% of forest area, which was not defined by any model (named not defined).  (Tables S1 and S2).
The total forest area occupied by the forest types in the Iberian Peninsula is 206.478 km 2 , which represents around 35% of the total area. Figure 4 presents the dominance of each forest type in the study area. Seven forest types represents around 54% of the plots, mainly dominated by pines, eucalypts and perennial oaks, with the following composition: Open forests of pines, eucalypts or other species (o.M.mix); Dense mediterranean shrubs with Pinus halepensis (S.250); Pinus sylvestris forests (T.Py); Quercus ilex forests (T.Az); Open Quercus ilex forests with low understory (o.M.Az); Other pines forests (T.Px); and Pinus pinaster forests with understory (M.Pb). On the other hand, the five forest types that have the lower area of occupation, with less than 1% of the area, are dominated by tree broadleaves deciduous species, conifers and invasive species. They are: Quercus robur forests (T.Qr); Closed riparian forests (c.T.Rpx); Other deciduous oaks forests (T.Qx); Other conifer forests (T.Rx); and Invasive species formations (T.Ax). 1 The codes have the following rules: the first letter is related with canopy cover, where c is closed forests, and o is open forests; the second letter is related to the composition, where S is shrubs dominance, M is mixed dominance, and T is tree dominance; the third letter is related to the dominant species (or shrub formation) and defined by the code of the species, or group of species (Tables S1 and S2).
The total forest area occupied by the forest types in the Iberian Peninsula is 206.478 km 2 , which represents around 35% of the total area. Figure 4 presents the dominance of each forest type in the study area. Seven forest types represents around 54% of the plots, mainly dominated by pines, eucalypts and perennial oaks, with the following composition: Open forests of pines, eucalypts or other species (o.M.mix); Dense mediterranean shrubs with Pinus halepensis (S.250); Pinus sylvestris forests (T.Py); Quercus ilex forests (T.Az); Open Quercus ilex forests with low understory (o.M.Az); Other pines forests (T.Px); and Pinus pinaster forests with understory (M.Pb). On the other hand, the five forest types that have the lower area of occupation, with less than 1% of the area, are dominated by tree broadleaves deciduous species, conifers and invasive species. They are: Quercus robur forests (T.Qr); Closed riparian forests (c.T.Rpx); Other deciduous oaks forests (T.Qx); Other conifer forests (T.Rx); and Invasive species formations (T.Ax).  For each forest type, the total coverage for each vertical layer was calculated and represented in the Figure 5. Layer seven (h7) includes litter cover percentage (marked in blue). The different layers indicate the dominant species (trees, species groups, or shrub formations) in terms of percentage of cover for each species, in a vertical space defined by seven layers. The result is a detailed description of almost all the study area supported by a complete visualization of the model's distribution. Specifically, the forest types dominated by tree species have more species coverage in the upper height classes (from 4 m over 16 m: h4 to h1), whereas those dominated by shrub formations have more species coverage in the lower height classes (h7 to h5).
The analysis of the forest types allows the identification of the most common species present in the Iberian Peninsula. The average of maximum vegetation cover values, computed for the vertical layers of each forest type, has been used as a parameter to compare the abundance of each species. The percentage of forest area characterized by each species or group of species has also been computed. A threshold of 3% of vegetation cover has been used for area of incidence and cover values to compute the average. indicate the dominant species (trees, species groups, or shrub formations) in terms of percentage of cover for each species, in a vertical space defined by seven layers. The result is a detailed description of almost all the study area supported by a complete visualization of the model's distribution. Specifically, the forest types dominated by tree species have more species coverage in the upper height classes (from 4 m over 16 m: h4 to h1), whereas those dominated by shrub formations have more species coverage in the lower height classes (h7 to h5). The analysis of the forest types allows the identification of the most common species present in the Iberian Peninsula. The average of maximum vegetation cover values, computed for the vertical layers of each forest type, has been used as a parameter to compare the abundance of each species. The percentage of forest area characterized by each species or group of species has also been computed. A threshold of 3% of vegetation cover has been used for area of incidence and cover values to compute the average. When we compared the forest types with the NFIs classification for Portugal and Spain based on the dominant species (Table 3), it is clear that the first group of forest types (oaks and broadleaves species) and the second group (pines and eucalypts) correspond well to NFI classification. The forest types in the third group, dominated by shrubs, are associated to different NFI classes.
In terms of species prevalence, shrub formations are generally more spread in the forest area of the Iberian Peninsula, however the abundance in terms of cover is lower in comparison with tree species. On the contrary, tree species are less spread over all the forest area but with higher abundance of cover. In fact, around 71% of the clusters are dominated by trees, mainly species from the genera Pinus and Quercus, and 21% are dominated by shrub formations with species of Ulex spp., Cytisus spp., and Cistus spp. The Q. ilex s.l. and P. pinaster are the common species of importance for both NFIs. These results are in accordance with the distributions species in the Iberian Peninsula, with larger areas of heathlands, oak and pine forests in the Atlantic region, sclerophyllous, semi-deciduous forests and shrublands in the Mediterranean region with different pine species in the mountainous areas [14,15]. The literature indicates that around 50% of the potential natural and semi-natural forest areas are now covered essentially by pine and/or eucalypt plantations [14,15]. In fact, our results indicate that one cluster is dominated by the mixture P. pinaster and Eucalyptus spp. (o.M.mix) and several others have in their composition pines and eucalypts mixed with other tree species and/or shrub formations, for example forest type c.S.240, or c.S.210 or even c.S.110 ( Figure 5). The final descriptions of each forest type were represented schematically in Figure 6.

Ecological Preferences of the Forest Types
The distribution of the forest types is, as expected, affected by latitude ( Figure S1) and elevation gradients ( Figure S2). Forest types dominated by shrub formations are located mainly in lower elevations but in different latitudinal gradients. Forest types with riparian species (c.T.Rx), C. sativa (c.T.Ct), F. sylvatica (c.T.Fa) and P. sylvestris (T.Py) are distributed mainly in the north of the Iberian The European Forest Types (EFTs) classification, developed by the European Environmental Agency [39], with some adjustments made later by Barbati et al. [28], established categories and types of ecologically distinct forest communities dominated by specific assemblages of trees. A comparison between our forest types classification and the main EFTs category and type levels were performed based  Table S8. The forest types established for the Iberian Peninsula are composed by forests and other wooded lands and do not have a direct equivalence with the EFTs. Eight out of fourteen EFTs categories seems to have correspondence with the Iberian Peninsula forest types dominated by trees. Forest types with high level of heterogeneity (e.g., Open Pinus sylvestris and other pines forests with understory; and Open forests of pines, eucalypts or other species) are too broad to be assigned to an EFT category. Other forest types that are not based on individual species, such as Other deciduous oaks forests and Other pines forests, were marked as not categorized. The forest types dominated by P. pinaster, with and without understory (M.Pb and T.Pb), have correspondence with two EFT category and type levels (2.7 Atlantic Maritime pine forest and 10.1 Mediterranean pine forest) due to both forest types growing under two different bioclimates (temperate Oceanic and Mediterranean pluviseasonal-Oceanic) in Portugal and Spain [10]. An EFT with category 8 (Thermophilous deciduous forest) is well represented in the Mediterranean region dominated mainly by Quercus species, with or without secondary trees of Acer, Fraxinus and Carpinus species. This category also includes Castanea sativa forest (type level 8.7), an important forest for fruit production, distributed mostly in the northern part of both countries [40]. The type level 8.3 that includes Quercus pyrenaica forests is mainly located in the Ibero-Atlantic zone of Portugal and Spain under humid and sub-humid bioclimate [41]. An EFT with category 9 (Broadleaved evergreen forest) is also an important forest type in the Mediterranean and Macaronesian regions dominated by sclerophlyllous or lauriphyllous trees, mainly Quercus species. This category is well represented in the Iberian Peninsula by several forest types dominated by evergreen broadleaves oaks ( Figure 6). Forest types dominated by Mediterranean pines (e.g., T.Py, o.M.Pm, c.M.Pa) have correspondence with the EFT category 10 (Coniferous forest of the Mediterranean, Anatolian and Macaronesian regions).

Ecological Preferences of the Forest Types
The distribution of the forest types is, as expected, affected by latitude ( Figure S1) and elevation gradients ( Figure S2 The results of the relation between forest types with wildfire probability were organized in three groups by increased percentage of burned plots within groups ( Figure S5). In the first group, the highest wildfire probability was found in the forest type with Mediterranean influence climate (M.Sb Quercus suber forests with understory). Forest types dominated by other oaks or other broadleaves species, with a low dominance of understory, had the least proportion burned. The second group present several forest types dominated by pines, and the ones dominated mainly by P. pinaster, and eucalypts with understory have the higher wildfire probability. In the third group, the forest type Dense Ulex shrubs with pine or eucalypts (c.S.240) have the highest wildfire probability of all forest types. This forest type has more than 50% of coverage in each of the lower canopy layers (h < 2 m) dominated by Mediterranean xerophyllous and gorse shrublands (species of the genera Erica, Calluna and Ulex).

Discussion
Information about the composition of the canopy layers of forests is still non-standardized in Europe and is often collected in a non-harmonized way across countries because of differences in NFIs definitions, sampling designs, plot configurations, measured variables and measurement protocols [4]. Differences in definitions are the most important issues when harmonising the inventory results and there is large variations in measurement thresholds [1,5,26]. Encompassing the area of the Iberian Peninsula and based on Portugal and Spain NFIs data, we proposed a methodology to harmonize and thus compare NFIs estimates across administrative borders [42]. The harmonized forest composition threshold of the canopy layers was applied to shrubs and other species that dominate the lower height classes (understory, classes h7, h6 and h5) and trees in the higher height classes. Similar thresholds are commonly used in studies of other biological groups [43]. To facilitate the harmonising process on stand composition, between countries, it is needed to establish functions that allow us to perform the comparison. The percentage of tree crown cover and height by species, normally present in NFIs, enables establishing these functions. However, other variables which are not that frequent, cover and average height of shrubs, crown length, crown diameter or litter thickness should be also recommended.
To facilitate the analysis, interpretation and reporting of forest data, we established a forest type classification to convert a large forest territory to more homogeneous ecologically units that could be directly related to ecological forest characteristics [27,28,44]. The result of the K-means analysis produced the average percentage coverage of species according to height class (pseudo-species) for each cluster, thus outlining the profiles of resulting forest types. The establishment of forest types based on the composition of the stands was already explored by other authors in Portugal [45,46]. These authors took the information from the Fourth Portuguese National Forest Inventory (PTNFI4), that presented, for the first time, the coverage of species according to height classes. These measurements were considered of major importance to perform an assessment of the Portuguese forest biodiversity and were included in the PTNFI4, after the third Ministerial Conference on the Protection of Forests in Europe that was held in Lisbon in 1998 [47]. A forest type classification that allows the description of the territory in smaller and more homogeneous ecological units is very useful for many applications, e.g., to support decision makers for sustainable forest management, to establish strategies of ecosystems and conservation goals on protected areas and the Natura 2000 network, to monitor the provision of ecosystem services and get information concerning biodiversity [27,28,44]. The harmonized forest types classification can solve the lack of information at spatial level, which is crucial to strengthening cross-sectorial and international cooperation in biodiversity conservation policies and initiatives [48]. Conservation programs at a large scale are challenging because of the different political, social and ecological conditions between countries and economic systems [49]. Spatially harmonized information with more detail at a national level is needed [26,27] and in this study the 28 forest types obtained are able to describe the composition of the forests in the Iberian Peninsula and their distribution in a clear and at the same time detailed way. The distribution maps of each forest type generally do not show sharp differences between the vegetation on the borders of Portugal and Spain. The bioclimatic map indicates several bioclimates that are common to both countries [10]. The analysis of the species composition and the distribution maps gave, in general, an unambiguous view of the characteristics of each forest type.
Forest type classifications are usually organized on the basis of the main dominant tree species [27] and NFIs classification in both Portugal and Spain follow that structure. From the NFI results, per country, the most important tree species in Portugal are P. pinaster, Eucalytus globulus Labill. and Quercus suber L. [50], whereas in Spain the most important tree species are Q. ilex s.l., Pinus halepensis Miller, Pinus sylvestris L., Quercus pyrenaica Willd. and P. pinaster [51]. However, our results indicate that, even if the most common forest types are dominated by trees (Table 3), as expected, shrub formations are also dominant in areas classified as dominated by trees in the Iberian Peninsula. The role of shrubs is often not considered but the productivity of understory vegetation is probably comparable to that of the trees and has important biodiversity considerations [24,26]. Shrubs and understory vegetation are important elements in ecosystems and can contribute to the recovery of deforested and degraded landscapes [25]. The harmonized data could also be used to estimate Mean Species Cover (MSC), a key forest indicator, with information based on shrub cover [2].
The European Forest Types (EFTs) have a finer level of division in term of tree species composition, the type level, that generally facilitate the cross-links with national forest types classification systems [28,44] and Natura 2000 forest habitats at the EU scale [39]. There are no complete equivalences between the forest types as defined in this work and the ETFs categories/types. Some approximate correspondences are suggested in the Table S8 but a full assessment is not possible with the existing information. The more common ETF categories corresponding to Iberian forest types were categories 8 (Thermophilous deciduous forest), 9 (Broadleaved evergreen forest) and 10 (Coniferous forests of the Mediterranean, Anatolian and Macaronesian regions). In Italy, a harmonized forest map based on dominant species found similar bridges between their forest types systems and the EFTs system [27]. Different forest ecosystems resources assessments could possibly converge on a common nomenclature system for Europe [25][26][27]49].
We have analyzed other plot features to better characterize forest types. Species composition, distribution and regeneration are affected by spatial and temporal factors like climate, topography, aspect, slope position and land use [52]. Latitude and elevation gradients influenced the distribution of the forest types in the Iberian Peninsula. For instance, Q. suber, a Mediterranean species with an understory composed of shrublands, needs an average temperature of around 15°C to thrive and cannot tolerate very low temperatures, which limits its northern and altitudinal range [53]. Other authors also found that those two drivers determine the structure and tree species composition in forests. For example, in Andean forests [54] tree stem density and basal area increases with elevation while species richness decreases. In the same study, stem density and species richness were found to decrease with latitude. In general, for the northern hemisphere, south-facing slopes receive more sunlight and become more xeric and warmer, supporting drought-resistant vegetation and are less conducive for tree growth, while north-facing slopes retain moisture and are cold and humid, supporting moisture-loving plants. Our results indicate that temperate species (like Fagus sylvatica) require a humid atmosphere with precipitation well distributed throughout the year and a well-drained soil [40], avoiding south and west facing slopes, whereas Atlantic species (like Pinus sp. and Ulex sp.) prefer coastal areas with west and northwest slope aspect [40].
In Mediterranean regions, fire is part of the ecological process [35] and an important element in ecosystems dynamics [36]. Particularly in southern Europe wildfire represent a major element in many forested ecosystems and one of the main important natural hazards [37]. Forest community and ecosystem properties are important for the development of a wildfire and fires are major drivers of understory vegetation composition [24]. Fuel models used to predict stands future conditions should employ variables such as stand density, species composition and vertical structure, values that are known with reasonable accuracy [55]. Studies regarding the characterization of the vegetation and their structure, and thus the fuel models, have been conducted in order to settle the fire models of a forest. For example, Fernandes [19] established the fuel models for mainland Portugal based on data from the PTNFI5. In our study, the wildfire probabilities throughout the forest types were also analyzed. The highest wildfire probability was found in forest type of Dense Ulex shrubs with pine and eucalypts (c.S.240). Mediterranean gorse shrublands are communities with high horizontal and vertical vegetation continuity, in which the proportion of fine dead fuel fractions with low moisture content is around 50% of the total present phytomass [56]. Other forest types composed by pines and eucalypts species (e.g., o.M.mix, T.Pb) also have a higher probability to burn. These tree species have been the most affected by wildfires in the Iberian Peninsula [57,58]. Reducing fire hazards, by decreasing fuel loads and their continuity, requires management of these fire-prone shrublands in order to increase ecosystem resilience. The future development of systems of equations for estimating fuel load by fractions (foliage, thin and thick branches, etc.) and physiological state (green or dead), for the different tree and shrub species, will allow the incorporation of vertical distribution of fuels to each vegetation group that is proposed in this study. This information would allow us to estimate the likelihood of occurrence of surface, passive crown or active crown fires in each formation for different environmental conditions by using fire behaviour simulators. The development of large scale maps with this characterization would be very useful for a first approximation in the optimization of fuel management resources.

Conclusions
In this study, we propose a methodology to harmonize the information provided by NFIs of Portugal and Spain on the composition of the canopy layers to make comparability possible. The harmonized data were used to classify the forest formations in the Iberian Peninsula. Using a clustering algorithm, it was possible to establish 28 forest types for Portugal and Spain and their geographical distribution. The final description and cartography represent a useful tool describing almost all the Iberian Peninsula's forest lands. The NFIs have proven to be a valuable source of information for the development of tools to aid decision-making and the assessment of forest composition and their relationship with ecological characteristics.
The results can be readily updated with datasets from more recent NFIs. The developed method can be easily applied to different European countries and forest types, especially in other Mediterranean countries.
Supplementary Materials: The following are available online at http://www.mdpi.com/1999-4907/11/11/1170/s1, Table S1. Common groups of tree species in the Portuguese NFI (PTNFI) and Spanish NFI (SNFI), Table S2. Common groups of shrub formations in the Portuguese NFI (PTNFI) and Spanish NFI (SNFI). Codes of the formations according to SFM25. For each formation there are the corresponding species in Portugal and Spain NFIs, Table S3. Original vertical structure used in the SNFI3 (h represents the height), Table S4. Parameter estimates and coefficient of determination for the species or groups of species obtained fitting Equation (2), Table S5. Parameter estimates and coefficient of determination for the species or groups of species obtained fitting Equation (4), Table S6. Characteristics of the sample trees used in the development of the crown length equations [59], Table S7. Characteristics of the sample trees used in the development of the crown diameter equations [59], Table S8. Correspondence among forest types (FTs) and European Forest Types (EFTs), Figure S1. Funding: This research was funded from the European Union's Horizon 2020 research and innovation programme, DIABOLO project "Distributed, integrated and harmonized forest information for bioeconomy outlooks", EU Grant Agreement no. 633464, from the EXTREME project "Influence of VOCs (volatile organic compounds) on the extreme behavior of forest fires" (PCIF/GFC/0078/2018), Foundation for Science and Technology, and from the VIS4FIRE project (RTA 2017-0042-C05-05), Spanish Ministry of Economy,Industry, and Competitiveness.