Forestry, Faculty of
Non UBC
DSpace
Forests 11 (4): 409 (2020)
Marin, Gheorghe
Strimbu, Vlad C.
Abrudan, Ioan V.
Strimbu, Bogdan M.
In many countries, National Forest Inventory (NFI) data is used to assess the variability of forest growth across the country. The identification of areas with similar growths provides the foundation for development of regional models. The objective of the present study is to identify areas with similar diameter and basal area growth using increment cores acquired by the NFI for the three main Romanian species: Norway spruce (Picea abies L. Karst), European beech (Fagus sylvatica L.), and Sessile oak (Quercus petraea (Matt.) Liebl.). We used 6536 increment cores with ages less than 100 years, a total of 427,635 rings. The country was divided in 21 non-overlapping ecoregions based on geomorphology, soil, geology and spatial contiguousness. Mixed models and multivariate analyses were used to assess the differences in annual dimeter at breast height and basal area growth among ecoregions. Irrespective of the species, the mixed models analysis revealed significant differences in growth between the ecoregions. However, some ecoregions were similar in terms of growth and could be aggregated. Multivariate analysis reinforced the difference between ecoregions and showed no temporal grouping for spruce and beech. Sessile oak growth was separated not only by ecoregions, but also by time, with some ecoregions being more prone to draught. Our study showed that countries of median size, such as Romania, could exhibit significant spatial differences in forest growth. Therefore, countrywide growth models incorporate too much variability to be considered operationally feasible. Furthermore, it is difficult to justify the current growth and yield models as a legal binding planning tool.
ArticleRegional Variability of the Romanian Main TreeSpecies Growth Using National Forest InventoryIncrement CoresGheorghe Marin 1,2, Vlad C. Strimbu 3 , Ioan V. Abrudan 1 and Bogdan M. Strimbu 4,*1 Faculty of Silviculture and Forest Engineering, Transilvania University Brasov, 500036 Brasov, Romania;ghmarin@roifn.ro (G.M.); abrudan@unitbv.ro (I.V.A.)2 National Institute for Research and Development in Silviculture, 077190 Voluntari, Romania3 Faculty of Forestry, University of British Columbia, Vancouver, BC V6T 1Z4, Canada;vlad.strimbu@alumni.ubc.ca4 College of Forestry, Oregon State University, Corvallis, OR 97331, USA* Correspondence: bogdan.strimbu@oregonstate.eduReceived: 24 February 2020; Accepted: 2 April 2020; Published: 6 April 2020Abstract: In many countries, National Forest Inventory (NFI) data is used to assess the variabilityof forest growth across the country. The identification of areas with similar growths provides thefoundation for development of regional models. The objective of the present study is to identify areaswith similar diameter and basal area growth using increment cores acquired by the NFI for the threemain Romanian species: Norway spruce (Picea abies L. Karst), European beech (Fagus sylvatica L.),and Sessile oak (Quercus petraea (Matt.) Liebl.). We used 6536 increment cores with ages less than 100years, a total of 427,635 rings. The country was divided in 21 non-overlapping ecoregions based ongeomorphology, soil, geology and spatial contiguousness. Mixed models and multivariate analyseswere used to assess the differences in annual dimeter at breast height and basal area growth amongecoregions. Irrespective of the species, the mixed models analysis revealed significant differencesin growth between the ecoregions. However, some ecoregions were similar in terms of growthand could be aggregated. Multivariate analysis reinforced the difference between ecoregions andshowed no temporal grouping for spruce and beech. Sessile oak growth was separated not only byecoregions, but also by time, with some ecoregions being more prone to draught. Our study showedthat countries of median size, such as Romania, could exhibit significant spatial differences in forestgrowth. Therefore, countrywide growth models incorporate too much variability to be consideredoperationally feasible. Furthermore, it is difficult to justify the current growth and yield models as alegal binding planning tool.Keywords: spatial variability; annual rings; National Forest Inventory; ecoregions; mixed models;Picea abies; Fagus sylvatica; Quercus petraea1. IntroductionThe development of national forest growth and yield models served the society needs of the lastcentury, but fails to address the current challenges posed by the rise in population, expected standardof life, and constant changes of the environment [1]. The perceived endless computational poweravailable today recommends the distillation of the national models, such that they are tailored tospecific questions or areas. Such detailed models address not only the society’s expectations fromthe forest, but also serve as the basis of sustainable forest management. Arguably, the next refiningmodeling level is the identification of areas for which distinct models can be developed. The naturalidentification of such areas are the ecoregions. A plethora of definitions are available for an ecoregion,Forests 2020, 11, 409; doi:10.3390/f11040409 www.mdpi.com/journal/forestsForests 2020, 11, 409 2 of 18most of them being based on the work of Bailey [2], which relies on the idea of a “unified climate,geology, geomorphology, soil, potential natural vegetation, and predominant land use” [3]. Eventhough the idea of an ecoregion is clear, their delineation is rather controversial as they are biasedtoward a particular characteristic or approach (i.e., qualitative or quantitative), do not necessarilyprovide a holistic perspective on ecosystems, and do not necessarily operate at a specific level [4].Therefore, in this study we refer to an ecoregion as an area “within which there is spatial coincidencein characteristics of geographical phenomena associated with differences in the quality, health, andintegrity of ecosystems”, according to Omernik [4]. Although more variables (especially climaticvariables) could be included in ecoregion delineation, we focused on the analysis of non-climatic onesas they provide a cleared boundary definition. Irrespective of the details, an ecoregion is a pivotallandscape element, which is mandatory in regional and national management.The presence of ecoregions in national planning is particularly important when significantlydifferent climates and geomorphological areas span the country. The size of the country is not necessarilyrelevant, as even smaller countries can exhibit a large number of ecoregions, with an impact in forestmanagement. Romania, for example, is almost evenly distributed among the main geomorphologicalareas (i.e., plains, hills, and mountains), and its vegetation is particularly diverse [5–7].A significant source of forest data is the national forest inventories (NFIs). Most of the NFIseither experienced significant changes in the last 25 years in terms of design, the case of the USA orFrance [8,9], or were established for the first time, the case of the Eastern European countries [10].The Romanian NFI is the only agency that supplies reliable national and regional level informationfor the United Nations Convention on Climate Change, Convention on Biological Diversity or theForest Resource Assessment of the United Nations Food and Agriculture Organization [11]. TheRomanian NFI is a continuous forest inventory, using a systematic sampling design [10]. The inventorycovers the entire country during a 5-year cycle and is conducted in two phases: the first phase isphoto-interpretation, which uses orthophoto maps on the scale 1:5000 taken by aerial flights, andthe second phase is the terrestrial measurements. The measurements in areas of permanent sampleplots are performed continuously every NFI cycle under multi-annual planning. The NFI dividesthe country into squares of 4 × 4 km. The 4 × 4 km squares are sub-divided into 16 squares each of1 × 1 km. In the 1 × 1 km square located in the southwest corner of each 4 × 4 km square there is apermanent sample cluster (PSC) situated. The PSC is square-shaped with a side of 250 m. The sides ofthe squares are north-south or east-west oriented. There are four circular sample plots (SPs) situatedin each corner of the PSC. The field data is collected in each circular SP. The density of the PSCs isdetermined according to the general landscape classification of mountains, hills, and plains. In themountain and hilly regions, the NFI grid has cells of 4 × 4 km, whereas in the plain region the cells are2 × 2 km (Figure 1). The differential design of Romanian NFI captures the relevant information whilealso considering the density of the forest within the general landscape. Altogether, the Romanian NFIhas 31,201 PSC and 124,804 SP.Forests 2020, 11, 409 3 of 18 Forests 2020, 11, x FOR PEER REVIEW 3 of 18 3 Figure 1. The PSC of the Romanian NFI overlaid on the general location of the forest, in green. The color of the PSC reflects the geomorphology: red for plains, blue for hills, and black for mountains. The spatial distribution of the PSCs does not reflect the actual density [12]. The Romanian NFI collected a large number of increment cores, which represent valuable information for the development, refinement and spatial delineation of regional growth models. Analyzing NFI growth data in combination with regional spatial units, such as the ecoregions, represents the next logical step in the development of improved growth and yield models. Therefore, the objective of this study is to identify the existence of spatial variation in tree growth for the three main forest species in Romania: Norway spruce (Picea abies (L) H. Karst), European beech (Fagus sylvatica L), and Sessile oak (Quercus petraea (Matt.) Liebl). Furthermore, if the presence of spatial variation is confirmed, different patterns of annual tree growth are to be examined based on the characteristics of the ecoregions. Even though other studies have used increment cores from inventory data to analyze forest growth [13–15], the use of NFI increment cores data on a large area for model development represents, to our knowledge, a novelty in the literature. 2. Methods 2.1. Romanian Ecoregions One of the main drivers that could be responsible for tree growth variation across the country, if differences in growth exist, are likely the ecoregions. In this study, an ecoregion is defined as a contiguous area with similar geology, soils and geomorphology. The information describing the drivers of an ecoregion was supplied by the Romanian Academy. Based on the four drivers, we divided the country into 21 ecoregions (Figure 2). The ecoregions vary in size from 161,904 ha, the Figure 1. The PSC of the Romanian NFI overlaid on the general location of the forest, in green. Thecolor of the PSC reflects the geomorphology: red for plains, blue for hills, and black for mountains. Thespatial distribution of the PSCs does not reflect the actual density [12].The Romanian NFI collected a large number of increment cores, which represent valuableinformation for the development, refinement and spatial delineation of regional growth models.Analyzing NFI growth data in combination with regional spatial units, such as the ecoregions,represents the next logical step in the development of improved growth and yield models. Therefore,the objective of this study is to identify the existence of spatial variation in tree growth for the threemain forest species in Romania: Norway spruce (Picea abies (L.) H. Karst), European beech (Fagussylvatica L.), and Sessile oak (Quercus petraea (Matt.) Liebl). Furthermore, if the presence of spatialvariation is confirmed, different patterns of annual tree growth are to be examined based on thecharacteristics of the ecoregions. Even though other studies have used increment cores from inventorydata to analyze forest growth [13–15], the use of NFI increment cores data on a large area for modeldevelopment represents, to our knowledge, a novelty in the literature.2. Methods2.1. Romanian EcoregionsOne of the main drivers that could be responsible for tree growth variation across the country,if differences in growth exist, are likely the ecoregions. In this study, an ecoregion is defined as acontiguous area with similar geology, soils and geomorphology. The information describing thedrivers of an ecoregion was supplied by the Romanian Academy. Based on the four drivers, wedivided the country into 21 ecoregions (Figure 2). The ecoregions vary in size from 161,904 ha, theMaramures Plateau (323), to 2.8 million ha, the Transylvania Plateau (324) (Table 1). Out of the 21Forests 2020, 11, 409 4 of 18countrywide ecoregions, only 15 were used in this study, as within the remaining six, the three specieswere not encountered.Table 1. Romanian Ecoregions for which the growth of Norway spruce, European beech or Sessile oakwas analyzed [16].EcoregionNameEcoregionCode Soil Geomorphology Geology Province Area [ha]MoldavianPlateau 121Chernozems, Phaezems,Luvisols, Fluvisols,GleysolsPlateauSands, Gravels,Marls, Loess,SandstonesMoldavia 1,352,075MoldavianHills 122Luvisols, EutricCambisols, Phaezems HillsMarls, CalcarousMarls, Sandstones,Conglomerates,TuffsMoldavia 806,488EasternCarpathians 131Eutric Cambisols,Dystric Cambisols,Entic podzols, Haplicpodzol, RendzicleptosolMountains Flysch deposits Moldavia 1,546,525Buzau-Vranceapiedmonts 221Regosols, Luvisols,Eutric Cambisols HillsSands, Clay-sands,Marls, Sandstones,LimestonesMoldavia 595,917Getic plateau 222Chernozems, Phaezems,Luvisols, EutricCambisols, Vertisols,Fluvisols, GleysolsPlateauMarls,Sandy-marls,Sandstones, Tuffs,Sands, Gravels,Loess,Muntenia 1,934,478Buzau-Vranceamountains 231Dystric CambisolsEutric Cambisols MountainsFlysch depositsConglomerates Muntenia 550,591East SouthernCarpathians 232Dystric Cambisols Enticpodzols, Haplic podzol,Rendzic leptosolMountainsSchists, Phyllites,Gneiss,Conglomerates,Limestones,Muntenia 1,028,094West SouthernCarpathians 233Dystric CambisolsEutric Cambisols, Enticpodzols, Haplic podzol,Rendzic leptosolMountainsLimestones,Schists, Phyllites,GranitesMuntenia 449,804Caras Hills 321 Luvisols, Vertisols HillsArgillaceousmarls, Sands,GravelsBanat 444,045Cris Hills 322 Luvisols, EutricCambisols HillsSands, Gravels,Loess deposits Transylvania 679,685Maramuresplateau 323Eutric Cambisols,Luvisols PlateauMarls, Tuffs,Sandstones Transylvania 161,904TransylvaniaPlateau 324Luvisols, Phaezems,Chernozems PlateauMarls, Calcarousmarls tuffs,Sandstones,ConglomeratesTransylvania 2,798,764BanatMountains 331Eutric Cambisols,Luvisols, DystricCambisols, Enticpodzols, RendzicleptosolMountainsSchists, Granites,Gabbro,LimestonesBanat 708,961WesternCarpathians 332Eutric Cambisols,Luvisols, DystricCambisols, RendzicleptosolMountainsSchists, Gneiss,Granites,Limestones,SandstonesTransylvania 1,054,126Volcanic ridge 333 Andosols, Entic podzols Mountains Lahar deposits,Andesites Transylvania 950,029Forests 2020, 11, 409 5 of 18 Forests 2020, 11, x FOR PEER REVIEW 5 of 18 5 Figure 2. Romanian ecoregions delineated according to spatial contiguity, geology, soils, and geomorphology maps of the Romanian Academy [16]. 2.2. Increment Cores The species considered in this study, namely Norway spruce, European beech and Sessile oak, not only cover approximately 60% of the Romanian forests, but are also the main species in terms of economic, ecological, and social importance. To identify the presence of spatial variation in the growth of the tree species, we have used the increment cores of the Romanian NFI. The cores are extracted at each 5-year measurement at breast height (1.3 m) parallel with the contour line. The number of trees cored from each SP depends on the number of species: when only one species was present, then 3–4 cores were extracted, otherwise 2–3 cores/species were collected. For each cored tree, only one increment core was extracted and the cored trees are selected randomly, though are conditioned by specific constraints. The trees must belong to the dominant and co-dominant canopy class [17] and should not exhibit exterior damage. The cores are further mounted on a solid support, sanded, scanned and the width of each ring is measured using a graphical procedure. The main information extracted from the increment cores is the width of each ring and their succession. To ensure validity of the comparisons, we have included only the increment cores that were extracted from even-aged pure species stands. In addition, to ensure compatibility among various ecoregions, only stands with no thinning or no active management were considered. A more detailed description of the data can be found in Marin et al. [16]. The data used in this study consists of 6536 increment cores from 1655 PSC the last one being collected in 2011. The total number of rings is 427,635, out of which 136,904 rings are for Norway spruce, 241,240 rings for European beech, and 49,491 rings for Sessile oak. To focus the study on ages that are operational relevant, only rings with maximum age of 99 years were included. The number Figure 2. Romanian ecoregions delineated according to spatial contiguity, geology, soils, andgeomorphology maps of the Romanian Academy [16].2.2. Incre ent CoresThe species considered in this study, namely Norway spruce, European beech and Sessile oak,not only cover approximately 60% of the Romanian forests, but are also the main species in terms ofeconomic, ecological, and social importance. To identify the presence of spatial variation in the growthof the tree species, we have used the increment cores of the Romanian NFI. The cores are extractedat each 5-year measurement at breast height (1.3 m) parallel with the contour line. The number oftrees cored from each SP depends on the number of species: when only one species was present, then3–4 cores were extracted, otherwise 2–3 cores/species were collected. For each cored tree, only oneincrement core was extracted and the cored trees are selected randomly, though are conditioned byspecific constraints. The trees must belong to the dominant and co-dominant canopy class [17] andshould not exhibit exterior damage. The cores are further mounted on a solid support, sanded, scannedand the width of each ring is measured using a graphical procedure. The main information extractedfrom the increment cores is the width of each ring and their succession. To ensure validity of thecomparisons, we have included only the increment cores that were extracted from even-aged purespecies stands. In addition, to ensure compatibility among various ecoregions, only stands with nothinning or no active management were considered. A more detailed description of the data can befound in Marin et al. [16].The data used in this study consists of 6536 increment cores from 1655 PSC the last one beingcollected in 2011. The total number of rings is 427,635, out of which 136,904 rings are for Norwayspruce, 241,240 rings for European beech, and 49,491 rings for Sessile oak. To focus the study onages that are operational relevant, only rings with maximum age of 99 years were included. Thenumber of rings (cores) for each ecoregion varies according to the species (Table 2), from 272 rings(4 cores) (i.e., European beech in the Maramures Plateau (323)) to 49,576 (986 cores) (i.e., Norwayspruce in the Eastern Carpathians (131)). Each ring was described with five attributes: ecoregion,species, year, age, and width. The year is the calendar year of each ring within the increment cores, ageForests 2020, 11, 409 6 of 18represents the age of each ring, and ring width (RW) is the ring width in micrometers. Diameter atbreast height (DBH) and basal area (BA) increments are derived from ring width, and are consideredmore appropriate to characterize tree growth than RW. To avoid values for DBH growth determinedby special environmental conditions, we have used only cores from trees located inside the naturaldistribution range. Therefore, for Norway spruce, we have only included increment cores fromthe mountainous regions in the analysis, for European beech only those from the hill, plateau andmountainous regions, and for Sessile oak only those from the hill and plateau regions.Table 2. Size and number of annual rings and increment cores across ecoregions and species.Number of Rings Number of Increment CoresEcoregion (code) NorwaySpruceEuropeanBeechSessileOakNorwaySpruceEuropeanBeechSessileOakMoldavian Plateau (121) 2881 6575 41 93Moldavian Hills (122) 8687 2385 112 29Eastern Carpathians (131) 49,576 2286 986 332Buzau-Vranceapiedmonts (221) 11,572 4400 142 87Getic plateau (222) 15,496 13,962 223 219Buzau-Vranceamountains (231) 15,749 24,447 336 292East SouthernCarpathians (232) 27,969 30,700 506 387West SouthernCarpathians (233) 4259 20,280 84 258Caras Hills (321) 1577 1419 20 17Cris Hills (322) 2379 2969 31 44Maramures plateau (323) 272 706 4 11Transylvania Plateau (324) 13,763 17,075 187 253Banat Mountains (331) 1654 36,007 47 427Western Carpathians (332) 9102 21,448 228 285Volcanic ridge (333) 28,595 23,445 568 287Total 136,904 241,240 49,491 2755 3028 7532.3. AnalysisTo assess the differences in annual growth among ecoregions, we have used mixed model andmultivariate procedures. Each type of procedure addresses a different facet of the relationship betweenecoregions and growth: the mixed model reflects the annual growth within a neighborhood of one year,whereas the multivariate considers the annual growth within the context of all years. In the absence ofmeteorological data, we implicitly considered the change in climate by constraining the analysis to thelast 50 calendar years (i.e., since 1960). We have chosen as threshold the year 1960, as multiple studiespointed a shift in climate occurring at the end of the six decade of the twentieth century, not only inEurope [18–20] but also in the world [21–23]. Climate is an additional major driver of tree growth,besides site, but in this study we focus only on site, which forces the usage of a narrow time intervalfor assessment of spatial variability.2.3.1. Mixed Models AnalysisTo assess the impact of the ecoregion on the radial growth at breast height we have used alinear mixed modelling approach covariates age and DBH for DBH increment, or age and BA for BAForests 2020, 11, 409 7 of 18increment (Equations (1) and (2)). The ecoregions were treated as fixed effects whereas the age wasconsidered as random effect, as the age or each ring depends on the year.DBH incrementyear = Ecoregion+DBHyear−1 +Ageyear (1)BA incrementyear = Ecoregion+ BAyear−1 +Ageyear (2)We introduced age in the model to account for difference in radial growth with age [24–26]. Weselected DBH and BA as covariate to account for the impact of the size of the tree on the growth.Given the temporal dependency of the measurements and the nested structure of the data, weimplemented a repeated measurements framework within the liner mixed model. Ecoregion andthe covariate were considered as a fixed effect, whereas age was considered random. Among thepossible covariance structures, we have considered the compound symmetry, variance components,and autoregressive of order 1 [27,28]. To identify the appropriate structure we have used the Akaike’sInformation Criterion (AIC) [29]. To formally determine the differences among ecoregions, if any, wehave used the estimated marginal means, as recommended by Gotway and Stroup [30]. To detectthe existence of residual autocorrelation, we used the Durbin–Watson test [31]. All analyses wereperformed in SAS 9.4 [32].2.3.2. Multivariate AnalysisTo simultaneously examine the spatial growth across all ecoregions we have used multivariateanalysis having as variables width of the annual rings. Therefore, the multivariate analysis should havefor each individual tree, and consequently each ecoregion, a series of at least 10 variables, dependingon its age. However, trees grow differently through time, which means that rings with the same agemust also be grouped by year. The requirement of having rings with the same age for each ecoregionand year since 1960 was not fulfilled. Therefore, instead of using individual trees, we have used groupsof trees that have the same age and year. Even with this consolidation of data, old rings will be rare,meaning that only few trees would be considered in analysis. Consequently, the ring age should bereduced to a value that ensures at least one tree per calendar year.Among the multivariate methods available for the analysis of the ring widths, we have usedprincipal component analysis (PCA), hierarchical cluster analysis, and canonical discrimination analysis.The PCA represents an orthogonal transformation of the variable matrix that is used for variablereduction and data exploration [32]. We performed PCA to investigate if there is any evidence ofgrouping the growth by ecoregion or if this grouping is determined by a subset of ring widths.Hierarchical cluster analysis has a similar objective as PCA because it hierarchically clusters theobservations in the dataset, but using an unsupervised method as the number of clusters is notestablished apriori. We have created clusters using the Ward method [33]. The cubic clustering criterionand the scree plot analysis of eigenvalues were employed in determining the number of clusters [34,35].We complemented the hierarchical cluster analysis with canonical discrimination analysis, whichgrouped the ecoregions starting with an a priori number of groups [36]. We selected the number ofgroups of ecoregions as the number of clusters identified by the cluster analysis. The difference amongthe ecoregions was tested with Wilks’ Lambda, Pillai’s Trace, the Hotelling–Lawley Trace and Roy’sGreatest Root [33]. We executed all the analyses in SAS 9.4 [32].3. Results3.1. Mixed Models AnalysisThe tree growth charts revealed an obvious impact of the ecoregion on the diameter incrementand basal area increment (Figure 3). Depending on the species and the variable measuring growth,there are ecoregions that clearly have different growth increments than the rest, such as the EasternCarpathians (131) for Norway spruce, the Maramures plateau (323) for European beech, or the CarasForests 2020, 11, 409 8 of 18hills (321) for Sessile oak for DBHi, or the Buzau-Vrancea mountains (231) for BAi. However, thereseem to be regions which do not exhibit different increments, such as the Moldavian hills (122), theBuzau-Vrancea piedmonts (221), and the Transylvania Plateau (324) for the DBHi of Sessile oak. Theinference on ecoregion impact on growth based on visual information was supported by the summarystatistics for each species and attribute, as the increments varies by more than 100% for DBHi and BAifor European beech and Sessile oak (Table 3). The descriptive statistics paint a compelling argument,considering that the variance does not necessarily vary among ecoregions (Table 3). The summarystatistics and the charts suggest different increments for some ecoregions, which support the formalanalysis of the differences among ecoregions with mixed models.Table 3. Estimates and standard errors (in parenthesis) for DBH increment (DBHi) and basal areaincrement (BAi) for Norway spruce, European beech and Sessile oak by ecoregion.Ecoregion NameNorway Spruce European Beech Sessile OakDBHi[mm] BAi [mm2]DBHi[mm] BAi [mm2]DBHi[mm] BAi [mm2]Moldavian Plateau (121) 2.62(1.655)86.20(136.97)3.06(1.282)87.46(82.572)Moldavian Hills (122) 6.65(1.419)223.12(168.848)7.7(1.696)343.7(158.178)Eastern Carpathians (131) 4.81(1.790)242.41(183.54)3.26(1.437)95.93(119.595)Buzau-Vranceapiedmonts (221)3.79(1.291)116.32(99.491)3.45(1.514)109.59(110.656)Getic plateau (222) 3.49(1.434)104.26(111.399)2.78(1.160)86.88(79.274)Buzau-Vranceamountains (231)5.21(1.796)259.18(174.2)2.68(1.266)79.63(93.75)East SouthernCarpathians (232)4.06(1.552)205.57(148.561)3.19(1.293)103.07(99.455)West SouthernCarpathians (233)3.87(1.663)195.75(162.814)3.21(1.214)93.28(89.638)Caras Hills (321) 1.23(0.939)35.37(66.591)1.08(0.197)51.65(22.047)Cris Hills (322) 2.63(1.414)123.64(129.655)4.15(1.498)144.01(116.801)Maramures plateau (323) 7.16(1.336)327.48(143.497)6.13(1.738)230.53(232.496)TransylvaniaPlateau (324)3.32(1.418)98.47(118.498)3.4(1.392)120.2(114.688)Banat Mountains (331) 5.52(1.848)218.39(172.912)2.64(1.331)82.91(97.022)WesternCarpathians (332)5.27(1.699)236.36(167.421)3.09(1.158)88.88(82.883)Volcanic ridge (333) 4.35(1.570)220.25(161.321)3.27(1.502)108.19(125.989)Forests 2020, 11, 409 9 of 18 Forests 2020, 11, x FOR PEER REVIEW 8 of 18 8 there are ecoregions that clearly have different growth increments than the rest, such as the Eastern Carpathians (131) for Norway spruce, the Maramures plateau (323) for European beech, or the Caras hills (321) for Sessile oak for DBHi, or the Buzau-Vrancea mountains (231) for BAi. However, there seem to be regions which do not exhibit different increments, such as the Moldavian hills (122), the Buzau-Vrancea piedmonts (221), and the Transylvania Plateau (324) for the DBHi of Sessile oak. The inference on ecoregion impact on growth based on visual information was supported by the summary statistics for each species and attribute, as the increments varies by more than 100% for DBHi and BAi for European beech and Sessile oak (Table 3). The descriptive statistics paint a compelling argument, considering that the variance does not necessarily vary among ecoregions (Table 3). The summary statistics and the charts suggest different increments for some ecoregions, which support the formal analysis of the differences among ecoregions with mixed models. Figure 3. DBH increment and basal area increment mean by ecoregion for Norway spruce (a,b), European beech (c,d) and Sessile oak (e,f) over time. Table 3. Estimates and standard errors (in parenthesis) for DBH increment (DBHi) and basal area increment (BAi) for Norway spruce, European beech and Sessile oak by ecoregion. Ecoregion Name Norway Spruce European Beech Sessile Oak DBHi$$ [mm] BAi [mm2] DBHi [mm] BAi [mm2] DBHi [mm] BAi [mm2] Moldavian Plateau (121) 2.62 (1.655) 86.20 (136.97) 3.06 (1.282) 87.46 (82.572) Moldavian Hills (122) 6.65 (1.419) 223.12$$ (168.848) 7.7 (1.696) 343.7 (158.178) Eastern Carpathians (131) 4.81$$ (1.790) 242.41 (183.54) 3.26 (1.437) 95.93 (119.595) Buzau-Vrancea piedmonts (221) 3.79$$ (1.291) 116.32 (99.491) 3.45 (1.514) 109.59 (110.656) Getic plateau (222) 3.49 (1.434) 104.26 (111.399) 2.78 (1.160) 86.88 (79.274) Buzau-Vrancea mountains (231) 5.21$$ (1.796) 259.18 (174.2) 2.68 (1.266) 79.63 (93.75) East Southern Carpathians (232) 4.06$$ (1.552) 205.57 (148.561) 3.19 (1.293) 103.07 (99.455) Figure 3. DBH increment and basal area increment mean by ecoregion for Norway spruce (a,b),European beech (c,d) and Sessile oak (e,f) over time.While visual inspection and the summary statistics supported the hypothesis of the study thatthere are differences in increments among ecoregions, its formal assessment could reveal the existenceof a multitude of different growth increments. Therefore, the inclusion of covariates in the analysisseemed justified, as the visual inspection of the relationships DBHi–DBH and BAi–BA revealed notonly differences among some ecoregions, but also a trend, particularly for BAi (Figure 4). Forests 2020, 11, x FOR PEER REVIEW 9 of 18 9 West Southern Carpathians (233) 3.87 (1.663) 195.75 (162.814) 3.21 (1.214) 93.28 (89.638) Caras Hills (321) 1.23 (0.939) 35.37 (66.591) 1.08 (0.197) 51.65 (22.047) Cris Hills (322) 2.63 (1.414) 123.64 (129.655) 4.15 (1.498) 144.01 (116.801) Maramures plateau (323) 7.16 (1.336) 327.48 (143.497) 6.13 (1.738) 230.53 (232.496) Transylvania Plateau (324) 3.32 (1.418) 98.47 (118.498) 3.4 (1.392) 120.2 (114.688) Banat Mountains (331) 5.52$$ (1.848) 218.39 (172.912) 2.64 (1.331) 82.91 (97.022) Western Carpathians (332) 5.27$$ (1.699) 236.36 (167.421) 3.09 (1.158) 88.88 (82.883) Volcanic ridge (333) 4.35$$ (1.570) 220.25 (161.321) 3.27 (1.502) 108.19 (125.989) While visual inspection and the summary statistics supported the hypothesis of the study that there are differences in increments among ecoregions, its formal assessment could reveal the existence of a multitude of different growth increments. Therefore, the inclusion of covariates in the analysis seemed justified, as the visual inspection of the relationships DBHi–DBH and BAi–BA revealed not only differences among some ecoregions, but also a trend, particularly for BAi (Figure 4). Figure 4. Mean DBH increment versus mean DBH and mean basal area increment versus mean basal area by ecoregion for Norway spruce (a,b), European beech (c,d) and Sessile oak (e,f). Irrespective of the species, the AIC identified the first order autoregressive structure as the most suitable covariance structure. For Norway spruce, the ecoregions were found to have a significant effect on annual DBH growth and basal area growth. The highest estimate for average DBH increment was found in the Banat Mountains (331) ecoregion, and the lowest in the West Southern Carpathians (233) ecoregion. The estimated marginal means showed that most of the ecoregions are significantly different from one another in terms of growth (Table 4), except for eight pairs. In total, we found evidence that the DBH growth of Norway spruce significantly differs in six distinct ecoregions: the Figure 4. Mean DBH incre ent versus ean DB and ean basal area increment versus mean basalarea by ecoregion for Norway spruce (a,b), European beech (c,d) and Sessile oak (e,f).Irrespective of the species, the I i tifi t first r er autoregressive structure as the mostsuitable covariance structure. For or ay s r ce, t e ecoregions were found to have a significanteffect on annual DBH growth and basal area gro th. The highest esti ate for average DBH incrementwas found in the Banat Mountains (331) ecoregion, and the lowest in the West Southern Carpathians(233) ecoregion. The estimated marginal means showed that most of the ecoregions are significantlydifferent from one another in terms of growth (Table 4), except for eight pairs. In total, we foundevidence that the DBH growth of Norway spruce significantly differs in six distinct ecoregions: theForests 2020, 11, 409 10 of 18Eastern Carpathians (131), Buzau-Vrancea mountains (231), East Southern Carpathians (232), WestSouthern Carpathians (233), Banat Mountains (331), and Volcanic Ridge (333). The seventh regionanalyzed, the Western Carpathians (332), was found to be similar with four other regions. There is noevidence that the results are artifacts [37], as the main assumption needed to be checked in repeatedmeasurements analysis, namely that residuals are white noise, was fulfilled (i.e., the Durbin–Watsontest indicated no autocorrelation).Table 4. Impact on growth of Norway Spruce by ecoregions assessed with the estimated marginalmeans. Ecoregions not significantly different are identified for DBH increment and basal area incrementwith DBHi or BAi, respectively. The empty cells indicate significant difference.Ecoregion 131 231 232 233 331 332 333131231 DBHi232 BAi233 DBHi/BAi DBHi DBHi331 DBHi DBHi332 DBHi333The growth of European beech also varies significantly by ecoregion, the highest average DBHincrement being found in the Maramures Plateau (323), and the lowest in the Caras Hills (321) (Table 3).However, not all regions were significantly different from one another, mirroring the Norway spruce.In the case of European beech, we noticed groups of ecoregions that have both similar DBH and BAincrements (Table 5). We found nine groups of such ecoregions: Moldavian Plateau (121) with EasternCarpathians (131); Eastern Carpathians (131) with Caras Hills (321) and Transylvanian Plateau (324);Getic Plateau (222) with East Southern Carpathians (232), Caras Hills (321) and Volcanic Ridge (333);West Southern Carpathians (233) with Western Carpathians (332); Caras Hills (321) with TransylvanianPlateau (324) and Volcanic Ridge (333). Another 11 ecoregion pairs were similar in terms of DBHincrement and four in terms of BA increment (Table 5). Similarly to Norway spruce, there was noevidence that the results are based on methods violating analytical assumptions, as the Durbin–Watsontest indicated no autocorrelation of the residuals.Sessile oak mirrored the findings for Norway spruce and European beech, the ecoregions havinga significant effect on tree growth (Table 6). The highest average growth was in the Moldavian Hills(122) ecoregion and the lowest in the Caras Hills (321) (Table 3). The regions exhibited mostly differentgrowth from one another, with the exception of three pairs: the Buzau-Vrancea piedmonts (221) andCaras Hills (321); Caras Hills (321) and Cris Hills (322); Caras Hills (321) and Transylvania Plateau (324)(Table 6). There are also other four other similar pairs, but only from a DBH increment perspective.The Durbin–Watson test also found no autocorrelation of the residuals. It can be concluded that thereare seven ecoregions spanning the Romanian areal of Sessile oak that significantly differ in tree growth:the Moldavian Plateau (121), Moldavian Hills (122), Buzau-Vrancea piedmonts (221), Getic Plateau(222), Cris Hills (322), Maramures plateau (323) and Transylvania Plateau (324). Caras Hills (321) wasfound similar in DBH or BA increment with four other ecoregions.Forests 2020, 11, 409 11 of 18Table 5. Impact on growth of European beech by ecoregions assessed with the estimated marginal means. Ecoregions not significantly different are identified for DBHincrement and basal area increment with DBHi or BAi, respectively. The empty cells indicate significant difference.Ecoregion 121 122 131 221 222 231 232 233 321 322 323 324 331 332 333121 DBHi/BAi DBHi122 DBHi DBHi131 DBHi/BAi DBHi DBHi/BAi221 BAi DBHi/BAi DBHi DBHi222 DBHi/BAi DBHi/BAi BAi DBHi/BAi231 BAi232 DBHi/BAi BAi DBHi233 DBHi/BAi321 DBHi/BAi DBHi DBHi/BAi322323324331 DBHi332333Forests 2020, 11, 409 12 of 18Table 6. Impact on growth of Sessile oak by ecoregions assessed with the estimated marginal means.Ecoregions not significantly different are identified for DBH increment and basal area increment withDBHi or BAi, respectively. The empty cells indicate significant difference.Ecoregion 121 122 221 222 321 322 323 324121122221 DBHi/BAi DBHi DBHi222 DBHi321 DBHi/BAi DBHi/BAi322 DBHi3233243.2. Multivariate AnalysisAs expected, a value less than 50 ensured at least one tree per calendar year and ecoregion.Depending on the species, the number of annual rings that was used in all multivariate analyses iseither 20 for sessile oak or 30 for Norway Spruce and European beech. Irrespective of the attributemeasuring the growth, DBHi or BAi, the three multivariate analyses supported the finding of mixedmodel investigation for all three species, as a significant separation among regions was revealed. ThePCA, for which the first two principal components explained at most 80% of the growth variation,which occurred for BAi of sessile oak, (i.e., 79% for Norway Spruce and 73% for European beech), didnot reveal an evident grouping according to the ecoregion (Figure 5). Therefore, globally, there areother aspects of tree growth which are more important in delineating growth than ecoregions. Forests 2020, 11, x FOR PEER REVIEW 12 of 18 12 Table 6. Impact on growth of Sessile oak by ecoregions assessed with the estimated marginal means. Ecoregions not significantly different are identified for DB incre ent and basal area increment with DBHi or BAi, respectively. The e pty cells indicate significant difference. Ecoregion 121 221 222 321 322 323 324 121 122 221 DBHi/BAi DBHi DBHi 222 DBHi 321 DBHi /BAi DBHi /BAi 322 DBHi 323 324 3.2. Multivariate Analysis As expected, a value less than 50 ensured at least one tree per calendar year and ecoregion. Depending on the species, the number of annual rings that was used in all multivariate analyses is either 20 for sessile oak or 30 for Norway Spruce and European beech. Irrespective of the attribute measuring the growth, DBHi or BAi, the three multivariate analyses supported the finding of mixed model investigation for all three species, as a significant separation among regions was revealed. The PCA, for which the first two principal components explained at most 80% of the growth variation, which occurred for BAi of sessile oak, (i.e., 79% for Norway Spruce and 73% for European beech), did not reveal an evident grouping according to the ecoregion (Figure 5). Therefore, globally, there are other aspects of tree growth which are more important in delineating growth than ecoregions. Figure 5. Grouping by ecoregions for Norway spruce (a), European beech (b), and Sessile oak (c) with principal components 1 and 2 as axes. Because similar results were obtained for DBHi and BAi, we have represented only the results for DBHi. According to the cubic clustering criterion, the hierarchical cluster analysis found three groups of ecoregions for all three species (Figure 6), regardless the attribute measuring growth. Mirroring PCA findings, the groups were not as clearly delineated as mixed models analysis, but the results show distinct growing patterns among the three areas as pointed by the minimum difference among the groups (Table 7). The clear separation of the ecoregion is the results of the amount of information covered by the first two or three eigenvectors, which are responsible for the groups, that is at least 78%. Furthermore, the difference among the groups was more than 15% (the case of Norway spruce), but almost 50% for European beech and Sessile oak (i.e., >45%). The hierarchical clustering revealed that Norway spruce growing in the ecoregions Eastern Carpathians (131), Western Carpathians (332) and the volcanic ridge (333) are similar and distinct from the rest. The results were less conclusive for the European beech and sessile oak, as a clear pattern did not emerge from the tree hierarchy. Figure 5. Grouping by ecoregions for Norway spruce (a), European beech (b), and Sessile oak (c) withprincipal components 1 and 2 as axes. Because similar results were obtained for DBHi and BAi, wehave represented only the results for DBHi.According to the cubic clustering criterion, the hierarchical cluster analysis found three groups ofecoregions for all three species (Figure 6), regardless the attribute measuring growth. Mirroring PCAfindings, the groups were not as clearly delineated as mixed models analysis, but the results showdistinct growing patterns among the three areas as pointed by the minimum difference among thegroups (Table 7). The clear separation of the ecoregion is the results of the amount of informationcovered by the first two or three eigenvectors, which are responsible for the groups, that is at least 78%.Furthermore, the difference among the groups was more than 15% (the case of Norway spruce), butalmost 50% for European beech and Sessile oak (i.e., >45%). The hierarchical clustering revealed thatNorway spruce growing in the ecoregions Eastern Carpathians (131), Western Carpathians (332) andthe volcanic ridge (333) are similar and distinct from the rest. The results were less conclusive for theEuropean beech and sessile oak, as a clear pattern did not emerge from the tree hierarchy.Forests 2020, 11, 409 13 of 18 Forests 2020, 11, x FOR PEER REVIEW 13 of 18 13 Figure 6. Identification of the number of groups of ecoregions using Ward’s Minimum Variance and Cubic Clustering Criterion for Norway spruce (a), European beech (b) and Sessile oak (c). Because similar results were obtained for DBHi or BAi, we have represented only the results for DBHi. Table 7. Grouping of ecoregions with similar growth according to the hierarchical cluster analysis. Species Groups of Ecoregions Minimum Difference among Groups % of Variation Explained DBHi BAi Norway spruce 3 0.29 83.7 91.8 European beech 3 0.25 70.1 83.2 Sessile oak 3 0.25 71.5 87.5 Among all the multivariate methods, the canonical discrimination analysis was the only one that showed an obvious difference in growth between ecoregions (Figure 7), regardless of the attribute (i.e., DBHi or BAi). All four multivariate tests indicated that for all three species, the class means vectors are not equal. It was clear that Norway spruce from the Western Carpathians (131) grow differently than the West Southern Carpathians (233), and both from the rest of the ecoregions (Figure 7a). Whereas PCA and hierarchical cluster analysis failed to identify ecoregions with distinct growth from European beech, the canonical discriminant analysis succeeded, as the Moldavian plateau (121) and Moldavian hills (122) exhibit different growing patterns from the rest and between them (Figure 7b). We also found a clear delineation of ecoregions by growth for sessile oak, with the Buzau-Vrancea Piedmonts (221) showing different DBHi and BAi than the Transylvania plateau (324), which are both distinct from the other ecoregions. It should be noticed that fewer ecoregions can be considered by the canonical discriminant analysis (Figure 7c), as the ecoregions with less trees, consequently less factorial combinations, are not included (e.g., in the case of sessile oak, which has only four ecoregions, compared with eight in the mixed model analysis). Figure 7. Grouping by ecoregions for Norway spruce (a), European beech (b), and Sessile oak (c) with canonical variables as axes. Because the same conclusion is reached for DBHi and BAi only the plots for DBHi are presented. Figure 6. Identification of the number of groups of ecoregions using Ward’s Minimum Variance andCubic Clustering Criterion for Norway spruce (a), European beech (b) and Sessile oak (c). Becausesimilar results were obtained for DBHi or BAi, we have represented only the results for DBHi.rouping of ecoregions with s mila growth according to the hierarchical cluster analysis.Species Groups ofEcoregionsMini um Differenceamong Groups% of Variation ExplainedDBHi BAiNorway spruce 3 0.29 83.7 91.8European beech 3 0.25 70.1 83.2Sessile oak 3 0.25 71.5 87.5Among all the multivariate methods, the canonical discrimination analysis was the only one thatshowed an obvious difference in growth between ecoregions (Figure 7), regardless of the attribute (i.e.,DBHi or BAi). All four multivariate tests indicated that for all three species, the class means vectors arenot equal. It was clear that Norway spruce from the Western Carpathians (131) grow differently thanthe West Southern Carpathians (233), and both from the rest of the eco gi ns (Figu 7a). Whereas PCAand hiera chical cluster analysis failed o identify ecoregi ns with distinct gro th from European beech,the canonical discriminant alysis succeeded, as the Moldavian plateau (121) and Moldavian hills (122)exhibit different growing patterns from the rest and between them (Figure 7b). We also found a cleardelineation of ecoregions by growth for sessile oak, with the Buzau-Vrancea Piedmonts (221) showingdifferent DBHi and BAi than the Transylvania plateau (324), which are both distinct from the otherecoregions. It should be noticed that fewer ecoregions can be considered by the canonical discriminantanalysis (Figure 7c), as the ecoregions with less trees, consequently less factorial combinations, are notincluded (e.g., in the case of sessile oak, which has only four ecoregions, compared with eight in themixed model analysis). Forests 2020, 11, x FOR PEER REVIEW 13 of 18 13 Figure 6. Identification of the nu c regions using Ward’s Minimum Variance and Cubic Clustering Criterion for Nor a s r c ( ), r eech (b) and Sessile oak (c). Because similar results were obtained for DBHi or B i, e a e re rese te o ly the results for DBHi. Table 7. Grouping of ecoregions with similar growth according to the hierarchical cluster analysis. Species Groups of Ecoregions Minimum Difference among Groups % of Variation Explained DBHi BAi Norway spruce 3 0.29 83.7 91.8 European beech 3 0.25 70.1 83.2 Sessile oak 3 0.25 71.5 87.5 Among all the multivariate methods, the canonical discrimination analysis was the only one that showed an obvious difference in growth between ecoregions (Figure 7), regardless of the attribute (i.e., DBHi or BAi). All four multivariate tests indicated that for all three species, the class means vectors are not equal. It was clear that Norway spruce from the Western Carpathians (131) grow differently than the West Southern Carpathians (233), and both from the rest of the ecoregions (Figure 7a). Whereas PCA and hierarchical cluster analysis failed to identify ecoregions with distinct growth from European beech, the canonical discriminant analysis succeeded, as the Moldavian plateau (121) and Moldavian hills (122) exhibit different growing patterns from the rest and between them (Figure 7b). We also found a clear delineation of ecoregions by growth for sessile oak, with the Buzau-Vrancea Piedmonts (221) showing differe t DBHi an BAi than the Transylvania plateau (324), which are both distinct from the other ecoregions. It should be noticed that fewer ecoregions can be consid red by the canonical discriminant analysis (Figure 7c), as the ecoregions wit less trees, consequ ntly less fa torial combinations, are not included (e.g., in the case of sessile oak, which has only four ecoregions, compared with eight in the mixed model analysis). Figure 7. Grouping by ecoregions for Norway spruce (a), European beech (b), and Sessile oak (c) with canonical variables as axes. Because the same conclusion is reached for DBHi and BAi only the plots for DBHi are presented. Figure 7. Grouping by ecoregions for Norway spruce (a), European beech (b), and Sessile oak (c) withcanonical variables as axes. Because the same conclusion is reached for DBHi and BAi only the plotsfor DBHi are presented.Forests 2020, 11, 409 14 of 184. DiscussionTree growth is influenced by a wide range of factors across a multitude of gradients [38]. Knowledgeabout this variation can lead to improved and localized forest growth models, with direct impacton forest management, forest operations and planning [39–41]. Our analysis was able to identifyspatial/regional variations in the annual growth of the main tree species even for a country close to themedian size (ranked 81 by IndexMundi), such as Romania.The mixed models analysis showed that growth significantly differs by ecoregion. For all threespecies, several ecoregions grouped, indicating that regionalization of various forest related attributesor models is required. In the case of Norway spruce, there was similarity in DBH and BA incrementbetween ecoregions Western Southern Carpathians (233) and Banat Mountains (331). One explanationfor the analogous growth can be geomorphology, as all regions were mountainous, have similar soilsand aspect (south facing) and are spatially adjacent. Spatial proximity can represent an important factorwhen comparing ecoregions because it can encompass numerous factors that influence growth. Amongthe most cited of those factors is climate, with effects at multiple levels in tree physiology [42–44],but others such as genetics, disturbances or localized stress events can also have an impact on treedevelopment [45]. The DBH growth in the Banat Mountains (331) was also not different from WesternCarpathians (332) and Volcanic ridge (333) ecoregions. This correspondence could be based ongeomorphological factors, as soils, geology and location are different or it can also be attributed toother confounding factors. DBH growth also did not differ between Western Southern Carpathians(233) and Western Carpathians (332) and Volcanic ridge (333). As the ecoregions Banat Mountains (331)and Volcanic ridge (333) show similar DBH increment and mostly BA increment, these two ecoregionsare good candidates for being joined from the perspective of Norway spruce growth. The highestaverage growth (Banat Mountains (331)) and the lowest (West Southern Carpathians (233)) occurredin a mountainous region with a similar soil composition, but a different altitude. Banat Mountains(331) are having a lower altitude with stepped leveling appearance and remnant limestone geology,whereas West Southern Carpathians (233) present higher altitudes, with leveled relief and limestonegeology. The altitude can represent a proxy for other variables influencing growth, such as temperature,precipitation or for areas prone to stress and disturbances.In the case of European beech, similar growth was found between pairs of plateau regions(ecoregions Buzau-Vrancea piedmonts (221), Getic plateau (222), Caras Hills (321), Transylvanianplateau (324)). As beech has its native range at these altitudes, it is expected that growth in theseecoregions is higher and also closer in terms of increments. However, the fact that these ecoregions arenot completely equivalent, illustrates that the spatial variability of growth is present even within samegeomorphology and that other factors can be involved in determining tree growth. These can be soilor geology, but also climatic variables or the presence of localized stress events. The average growthin the Maramures Plateau (323) is different from all the other ecoregions, with the highest basal areaincrement estimate. One possible explanation for this difference could be attributed to the combinedeffects of growing conditions and other factors such as regional climate or natural distribution of thegenetic material. The growth similarity between Eastern Carpathians (131) and some plateau regions(Moldavian Plateau (121), Caras Hills (321), Transylvania Plateau (324)) indicates again that there aremultiple factors influencing growth, and that geomorphology is not the only defining factor.For Sessile oak, only a few ecoregions showed growth similarities. The Caras Hills (321) ecoregionexhibited similar growth with four other ecoregions where Sessile oak is present (Buzau-Vranceapiedmonts (221), Getic plateau (222), Cris Hills (322) and Transylvania Plateau (324)), all being locatedat the foothills of the mountains. This indicates that this foothill ecoregion represents a location ofaverage growth when compared with other Sessile oak ecoregions, and could be combined with theothers. Highest average growth occurred in a hilly ecoregion (Moldavian Hills (122)) and a muchlower one in the neighboring one, the Moldavian Plateau (121). The unambiguous difference betweenthese two regions underlines the spatial differentiation in tree growth for Sessile oak. As in NorwayForests 2020, 11, 409 15 of 18spruce, the altitude represents a potential differentiation factor, with the plateau ecoregion beingaround 200–500 m high and the hills around 500–900 m high.Two models for tree growth were developed in this study, and ecoregions could be differentiatedby both models. Some authors [46,47] suggest that both DBH and BA increment are appropriate ingrowth modelling, while Russell [48] suggests that DBH increment might decrease the model errormore than BA. The results of this study showed that, for all three species and especially for Europeanbeech and Sessile oak, both the DBH increment and BA increment models converged towards sameresults in a high number of cases. This is a clear indication that ecoregions have an effect on treegrowth, irrespective of the type of variable used to measure increment. There were instances whereonly DBH increment indicated a growth differentiation between ecoregions, as this variable might bemore suitable for increment modelling [48].The multivariate analysis examined simultaneously all the DBHi or BAi, looking at the growthpattern along the life of each individual tree. The PCA and cluster analysis suggested that DBHicannot evidently group ecoregions and cannot be grouped by certain periods of time. This means thatwithin the same species, variation of growth in time did not represent a grouping factor. The canonicaldiscrimination analysis indicated that there is growth difference by ecoregion for both Norway spruceand European beech. However, for Sessile oak, the canonical variables have discriminatory powerfor classifying growth into ecoregions, irrespective of DBHi or BAi. As the canonical variables area linear combination of yearly growth, DBHi in certain years could determine the inclusion in anecoregion. Sessile oak has its native range in warm areas which are more prone to draught and havingwater as limiting factor. Although further analysis in combination with climate data is required, it ispossible that certain ecoregions were more influenced by changes in climatic conditions, which wereinstrumental in this growth differentiation of ecoregions.One main finding of the mixed models analysis was that geomorphology plays a major role in thespatial variation of growth. Ecoregions were found to have similar growth increments, even whensoil, geology or location differ. The geomorphology was found to be relevant for all three species.One explanation for this finding is that geomorphology changes the specific soil characteristics [49]and indirectly influences growth. It is also possible that geomorphology has driven the historicalhuman impact on forest cover and growth (i.e., higher impact in flatter areas and lower in mountainsones), and therefore influencing in this alternate manner the vegetation development [50]. Anotherexplanation can be that climate conditions present in that geomorphology type have a strong impacton tree growth [44,51,52], which can be a stronger driver than soil [53]. Non-climatic variables such asgeomorphology can represent a proxy for different climatic characteristics [54,55], as spatially adjacentecoregions with different altitudes exhibited very different growth increments. However, even whenconsidering only non-climate variables, spatial variability was evident.5. ConclusionsThe spatial variation of growth for tree species represents an essential aspect of forest managementplanning, modelling and reporting [56]. This study aimed at determining the existence of variation inDBH and BA growth over the area of a medium sized country located in a temperate region. A seriesof 21 ecoregions summarized the ecological conditions for the entire country.Using mixed models and multivariate analysis on 6536 increment cores from the RomanianNFI, we were able to establish that tree growth is significantly different between various ecoregionsfor Norway spruce, European beech and Sessile oak. Numerous factors such as geology, soil andgeomorphology are responsible for the difference in growth, with the latter being frequently present inthe differentiation of ecoregions. An important finding of this study was that global, country wide,growth models incorporate far too much variability from an operational perspective. They can beused for country-wide resource estimations, but it is hard to justify their usage as a legal bindingplanning tool.Forests 2020, 11, 409 16 of 18As the change in climate also influences the vegetation growth, and consequently, the type andamount of products and services supplied by each ecosystem, further dendrological studies thatinclude climate and temporal analysis are needed for a better understanding of tree growth variability.Author Contributions: G.M. organized data collection and framed the study, V.C.S. interpreted and wrote themanuscript, I.V.A. wrote portions of the manuscript, B.M.S. analyzed, interpreted, and wrote portions of themanuscript. All authors have read and agreed to the published version of the manuscript.Acknowledgments: This research was partially funded by the National Institute of Food and Agriculture, U.S.Department of Agriculture, grant number 2019-67019-29462, the McIntire Stennis project OREZ-FERM-875, andthe Romanian National Authority for Scientific Research and Innovation, CCCDI—UEFISCDI, project numberERANET-FACCE ERAGAS—FORCLIMIT (82/2017), within PNCDI II.Conflicts of Interest: The authors declare no conflict of interest.References1. Zianis, D.; Muukkonen, P.; Mäkipää, R.; Mencuccini, M. Biomass and Stem Volume Equations for Tree Species inEurope; Silva Fennica Monographs 4: Tampere, Finland, 2005.2. Bailey, R.G. Delineation of ecosystem regions. Environ. Manag. 1983, 7, 365–373. [CrossRef]3. Hughes, R.M.; Omernik, J.M. Ecological regions (ecoregions). In Environmental Geology; Springer: Dordrecht,The Netherlands, 1999; pp. 155–159. ISBN 978-1-4020-4494-6.4. Omernik, J.M. Perspectives on the nature and definition of ecological regions. Environ. Manag. 2004, 34,S27–S38. [CrossRef] [PubMed]5. Indreica, A.; Turtureanu, P.D.; Szabó, A.; Irimia, I. Romanian forest database: A phytosociological archive ofwoody vegetation. Phytocoenologia 2017, 47, 389–392. [CrossRef]6. Šijacˇic´-Nikolic´, M.; Milovanovic´, J.; Nonic, M. (Eds.) Forests of Southeast Europe under a Changing Climate:Conservation of Genetic Resources; Advances in Global Change Research; Springer International Publishing:Berlin/Heidelberg, Germany, 2019; ISBN 978-3-319-95266-6.7. Stapf, O. Conspectul Florei României. Nature 1899, 59, 221. [CrossRef]8. Bechtold, W.A.; Patterson, P.L.The Enhanced Forest Inventory andAnalysis Program -National SamplingDesign andEstimation Procedures; U.S. Department of Agriculture, Forest Service, Southern Research Station: Asheville,NC, USA, 2005; p. 85.9. Vidal, C.; Belouard, T.; Herve, J.C.; Rober, N.; Wolsak, J. A New Flexible Forest Inventory in France; USDAForest Service: Portland, ME, USA, 2007; pp. 67–73.10. Tomppo, E.; Gschwantner, T.; Lawrence, M.; McRoberts, R.E. (Eds.) National Forest Inventories: Pathways forCommon Reporting; Springer: Dordrecht, The Netherlands, 2010; ISBN 978-90-481-3232-4.11. IFN National Forest Inventory. Available online: https://roifn.ro/site/en/ (accessed on 6 August 2019).12. Marin, G.; Bouriaud, O.; Dumitru, M.; Nitu, D. Romania. In National Forest Inventories: Pathways for CommonReporting; Springer: Heidelberg, Germany, 2010; pp. 473–480.13. Charru, M.; Seynave, I.; Morneau, F.; Bontemps, J.-D. Recent changes in forest productivity: An analysisof national forest inventory data for common beech (Fagus sylvatica L.) in north-eastern France. For. Ecol.Manag. 2010, 260, 864–874. [CrossRef]14. Evans, M.E.K.; Falk, D.A.; Arizpe, A.; Swetnam, T.L.; Babst, F.; Holsinger, K.E. Fusing tree-ring and forestinventory data to infer influences on tree growth. Ecosphere 2017, 8, e01889. [CrossRef]15. Metsaranta, J.M. Dendrochronological procedures improve the precision and accuracy of tree and stand ageestimates in the western Canadian boreal forest. For. Ecol. Manag. 2020, 457, 117657. [CrossRef]16. Marin, G.; Abrudan, I.; Strimbu, B.M. Increment cores of the National Forest Inventory from Romania. Math.Comput. For. Nat. Resour. Sci. 2019, 11, 294.17. Nyland, R.D. Silvicuture: Concepts and Applications; McGraw-Hill: New York, NY, USA, 1996.18. Kuglitsch, F.G.; Toreti, A.; Xoplaki, E.; Della-Marta, P.M.; Zerefos, C.S.; Türkes¸, M.; Luterbacher, J. Heat wavechanges in the eastern Mediterranean since 1960. Geophys. Res. Lett. 2010, 37. [CrossRef]19. Zuidhoff, F.S.; Kolstrup, E. Changes in palsa distribution in relation to climate change in Laivadalen, northernSweden, especially 1960–1997. Permafr. Periglac. Process. 2000, 11, 55–69. [CrossRef]20. Gregory, P.J.; Marshall, B. Attribution of climate change: A methodology to estimate the potential contributionto increases in potato yield in Scotland since 1960. Glob. Chang. Biol. 2012, 18, 1372–1388. [CrossRef]Forests 2020, 11, 409 17 of 1821. Wang, S.; Zhang, M.; Li, Z.; Wang, F.; Li, H.; Li, Y.; Huang, X. Glacier area variation and climate change in theChinese Tianshan Mountains since 1960. J. Geogr. Sci. 2011, 21, 263–273. [CrossRef]22. Yue, T.X.; Fan, Z.M.; Liu, J.Y. Changes of major terrestrial ecosystems in China since 1960. Glob. Planet. Chang.2005, 48, 287–302. [CrossRef]23. Posmentier, E.S.; Cane, M.A.; Zebiak, S.E. Tropical pacific climate trends since 1960. J. Clim. 1989, 2, 731–736.[CrossRef]24. Husch, B.; Beers, T.W.; Kershaw, J.A. Forest Mensuration, 4th ed.; John Wiley and Sons: Hoboken, NJ,USA, 2002.25. Pretzsch, H. Forest Dynamics, Growth and Yield; Springer: Berlin, Germany, 2009.26. Weiskittel, A.R.; Hann, D.W.; Kershaw, J.A.; Vanclay, J.K. Forest Growth and Yield Modeling; Wiley-Blackwell:Chichester, UK, 2011; Volume 430.27. Crowder, M.J.; Hand, D.J. Analysis of Repeated Measures; Monographs on Statistics and Applied Probability;Chapman and Hall: London, UK, 1990.28. Fitzmaurice, G.M.; Laird, N.M.; Ware, J.H. Applied Longitudinal Analysis; Wiley: Hoboken, NJ, USA, 2004;ISBN 978-0-471-21487-8.29. Pinheiro, J.; Bates, D. Mixed Effects Models in S and S-PLUS; Springer: New York, NY, USA, 2000.30. Gotway, C.A.; Stroup, W.W. A generalized linear model approach to spatial data analysis and prediction. J.Agric. Biol. Environ. Stat. 1997, 2, 157–178. [CrossRef]31. Neter, J.; Kutner, M.H.; Nachtsheim, C.J.; Wasserman, W. Applied Linear Statistical Models; WCB McGraw-Hill:Boston, MA, USA, 1996.32. SAS Institute. SAS®9.4; SAS Institute Inc.: Cary, NC, USA, 2013.33. Rencher, A.C. Methods of Multivariate Analysis; John Wiley and Sons: New York, NY, USA, 2002.34. Hardle, W.; Simar, L. Applied Multivariate Statistical Analysis; Springer: New York, NY, USA, 2003.35. Tabachnick, B.G.; Fidell, L.S.UsingMultivariate Statistics; Allyn and Bacon: Needham Heights, MA, USA, 2001.36. Strimbu, B.M.; Hickey, G.M.; Strimbu, V.G.; Innes, J.L. On the use of statistical tests with non-normallydistributed data in landscape change detection. For. Sci. 2009, 55, 72–83.37. Seppelt, R.; Richter, O. “It was an artefact not the result”: A note on systems dynamic model developmenttools. Environ. Model. Softw. 2005, 20, 1543–1548. [CrossRef]38. Oliver, C.D.; Larson, B.C. Forest Stand Dynamics, updated ed.; John Wiley & Sons, Inc.: New York, NY, USA,1996; ISBN 978-0-471-13833-4.39. Hasenauer, H. Concepts within tree growth modeling. In Sustainable Forest Management: Growth Modelsfor Europe; Hasenauer, H., Ed.; Springer Berlin Heidelberg: Berlin/Heidelberg, Germany, 2006; pp. 3–17.ISBN 978-3-540-31304-5.40. Strimbu, V.C.; Bokalo, M.; Comeau, P.G. Deterministic models of growth and mortality for jack pine in borealforests of western Canada. Forests 2017, 8, 410. [CrossRef]41. Taylor, A.R.; Chen, H.Y.H.; VanDamme, L. A review of forest succession models and their suitability forforest management planning. For. Sci. 2009, 55, 23–36.42. Savva, Y.; Oleksyn, J.; Reich, P.B.; Tjoelker, M.G.; Vaganov, E.A.; Modrzynski, J. Interannual growth responseof Norway spruce to climate along an altitudinal gradient in the Tatra Mountains, Poland. Trees 2006, 20,735–746. [CrossRef]43. Andreassen, K.; Solberg, S.; Tveito, O.E.; Lystad, S.L. Regional differences in climatic responses of Norwayspruce (Picea abies L. Karst) growth in Norway. For. Ecol. Manag. 2006, 222, 211–221. [CrossRef]44. Kramer, K.; Leinonen, I.; Loustau, D. The importance of phenology for the evaluation of impact of climatechange on growth of boreal, temperate and Mediterranean forests ecosystems: An overview. Int. J. Biometeorol.2000, 44, 67–75. [CrossRef] [PubMed]45. Zhang, Q.; Alfaro, R.I.; Hebda, R.J. Dendroecological studies of tree growth, climate and spruce beetleoutbreaks in Central British Columbia, Canada. For. Ecol. Manag. 1999, 121, 215–225. [CrossRef]46. Shifley, S.R. A Generalized System of Models Forecasting Central States Tree Growth; Research Paper NC-279; USDepartment of Agriculture, Forest Service, North Central Forest Experiment Station: St. Paul, MN, USA,1987; Volume 279.47. West, P.W. Use of diameter increment and basal area increment in tree growth studies. Can. J. For. Res. 1980,10, 71–77. [CrossRef]Forests 2020, 11, 409 18 of 1848. Russell, M.B.; Weiskittel, A.R.; Kershaw, J.A., Jr. Comparing strategies for modeling individual-tree heightand height-to-crown base increment in mixed-species Acadian forests of northeastern North America. Eur. J.For. Res. 2014, 133, 1121–1135. [CrossRef]49. Griffiths, R.P.; Madritch, M.D.; Swanson, A.K. The effects of topography on forest soil characteristics in theOregon Cascade Mountains (USA): Implications for the effects of climate change on soil properties. For. Ecol.Manag. 2009, 257, 1–7. [CrossRef]50. Sandel, B.; Svenning, J.-C. Human impacts drive a global topographic signature in tree cover. Nat. Commun.2013, 4, 2474. [CrossRef]51. Laubhann, D.; Sterba, H.; Reinds, G.J.; De Vries, W. The impact of atmospheric deposition and climate onforest growth in European monitoring plots: An individual tree growth model. For. Ecol. Manag. 2009, 258,1751–1761. [CrossRef]52. Scharnweber, T.; Manthey, M.; Criegee, C.; Bauwe, A.; Schröder, C.; Wilmking, M. Drought matters–Decliningprecipitation influences growth of Fagus sylvatica L. and Quercus robur L. in north-eastern Germany. For. Ecol.Manag. 2011, 262, 947–961. [CrossRef]53. Toledo, M.; Poorter, L.; Peña-Claros, M.; Alarcón, A.; Balcázar, J.; Leaño, C.; Licona, J.C.; Llanque, O.;Vroomans, V.; Zuidema, P.; et al. Climate is a stronger driver of tree and forest growth rates than soil anddisturbance. J. Ecol. 2011, 99, 254–264. [CrossRef]54. Adams, H.R.; Barnard, H.R.; Loomis, A.K. Topography alters tree growth–climate relationships in a semi-aridforested catchment. Ecosphere 2014, 5, art148. [CrossRef]55. Tardif, J.; Camarero, J.J.; Ribas, M.; Gutiérrez, E. Spatiotemporal variability in tree growth in the centralpyrenees: climatic and site influences. Ecol. Monogr. 2003, 73, 241–257. [CrossRef]56. Spiecker, H. Tree rings and forest management in Europe. Dendrochronologia 2002, 20, 191–202. [CrossRef]© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open accessarticle distributed under the terms and conditions of the Creative Commons Attribution(CC BY) license (http://creativecommons.org/licenses/by/4.0/).
Article
eng
Reviewed
Vancouver : University of British Columbia Library
Multidisciplinary Digital Publishing Institute
CC BY 4.0
Faculty
spatial variability
annual rings
National Forest Inventory
ecoregions
mixed models
Picea abies
Fagus sylvatica
Quercus petraea
Regional Variability of the Romanian Main Tree Species Growth Using National Forest Inventory Increment Cores
Text