A model of shrub biomass accumulation as a tool to support management of Portuguese forests

Abstract

The evaluation of forest fuel loading is required by most fire management activities. However, the consideration of shrub biomass for forest planning purposes has been limited by the inability to predict its growth and accumulation. The main objective of this study was to model shrub biomass over time under a tree canopy to be able to include shrub management in fire risk mitigation plans. Data for this purpose was obtained from the 4^th and 5^th Portuguese National Forest Inventories. Five biologically realistic models were built to describe shrub biomass accumulation in Portuguese forests. The selected model indicates that maximum biomass is affected by stand basal area and the percentage of resprouting shrub species in the stand. Biomass growth rate was clearly affected by the regeneration strategies of shrubs in combination with climatic conditions (mean annual temperature). The model can be used in the accumulation form for initialization purposes or in one of the two alternative difference forms to project observed shrub biomass. The acquired ability to estimate shrub biomass facilitates its inclusion in forest growth models and simulators and will contribute to more accurate estimates of fire behaviour characteristics and stored carbon. This is instrumental to improve decision-making in forest management plans that integrate fire risk, namely to schedule understory fuel treatments.

Introduction

During the last decades, significant land use change took place in Portugal and elsewhere in the Mediterranean region. Many marginal agricultural or grazing lands were either abandoned or afforested. Natural succession led to changes in vegetation structure and composition where agricultural activities ceased, thus contributing to the expansion of shrubland, woodland and forests with a well-developed shrub understory ([30]). These changes resulted in higher carbon stocks as well as in more flammable ecosystems prone to large and high-severity fires ([60], [15]).

Fuel dynamics refers to the structural and temporal modifications undergone by a fuel layer or fuel complex. Shrub accumulation models could assist in forecasting the dynamics of biomass and carbon storage. Many modeling of fuel dynamics follows the simple model by Olson ([54]) that describes the relationship between production and decomposition as a modified exponential function that flattens out to a plateau. Other studies describe fuel and shrub dynamics by time-dependent models of forest fire hazard ([37]). However, shrub biomass accumulation information for Mediterranean areas is very limited. Few studies addressed the temporal dynamics of shrub structure and/or biomass in shrublands ([3]), which are expected to be different under a forest canopy, due to competition for resources (i.e., light, water). Hence, little attention has been given to understory vegetation, likely due to its limited economic importance. Nonetheless, the ecological significance of the understory is high, since it plays an important role on nutrient cycles, carbon storage and fire hazard.

Currently available carbon models still lack details on biomass dynamics, which in turn affect the calculation of these processes. A recent study by Rosa et al. ([67]) to estimate pyrogenic emissions of greenhouse gases, aerosols and other trace gases from wildfires in Portugal identified shrub biomass as the variable with the greatest impact on the uncertainty inherent in such estimates. Therefore, it is essential to improve the assessment of forest biomass, including its spatial and temporal variation.

In the Mediterranean region, fire is one of the most important factors affecting forest ecosystems, both ecologically and economically ([63]). Higher shrub loading implies higher flammability, likelihood of crowning fire, and difficulty in fire control ([69], [28]). Fernandes et al. ([26]) observed differences in fire behavior and severity among maritime pine (Pinus pinaster) plots depending on fuel age (i.e., time since last treatment) and the presence or absence of surface fuel treatments. Furthermore, recent research aimed at developing wildfire occurrence models in Portugal included the understory shrubs biomass as a significant variable. Indeed, shrubs have a large impact on fire risk with obvious implications to forest planning ([33], [34], [48], [8]).

The integration of wildfire risk in forest management planning depends on the continuously changing variables related with fuel dynamics (e.g., tree and shrub growth) and stand management ([36], [35]). Therefore, it is very important to obtain information as much accurate as possible on all key variables (e.g., shrub growth) affecting the likelihood and severity of fire over time. The lack of tools to project shrub growth over time has hindered its inclusion in the forest management planning ([35]).

A broad range of growth modeling techniques applied to forest ecosystems has been reported in the scientific literature. However, biomass growth and yield information of understory is scarce. The usefulness of shrub biomass models in forest planning depends on the input information they need and whether the future values of predictors can be estimated with reasonable accuracy. Moreover, an increasing amount of information is being collected in forest inventories focused on sustainability and biodiversity. For example, the Portuguese National Forest Inventory (NFI) systematically collects plot-level information on the shrubs type (i.e., species and regeneration mechanism), their ground cover and mean height. This information may be used to develop a shrub biomass build-up model, which would contribute to more accurate estimates of fire behavior and may quantify the impact of silvicultural treatments on the probability of wildfire occurrence. It can also improve decision-making in forest management, especially taking into account the risk of forest fires ([31], [32]).

The aim of this study is to develop a model to describe the temporal dynamics of shrub biomass in the forest understory of Portuguese forests. One hypothesis is that maximum potential biomass, defined by the model asymptote, is affected by stand variables (e.g., basal area) and shrub regeneration types (resprouter vs. non-resprouter). Considering the biomass growth rate, our hypothesis is that it is affected by stand and shrub characteristics, as well as by site conditions and climate. This model may then be used to predict fuel accumulation as well as to update carbon stock inventories and may also be instrumental to include that information in forest management systems aiming at reducing the risk of fire.

Material and Methods

Study area

Mainland Portugal (Fig. 1) is located in southwester Europe at latitudes of 37° N to 42° N and longitudes of 6° W to 10° W. Climate is warm and dry during summer, and cool and wet in winter ([40]). Mean annual temperature and precipitation follow a gradient of increasing temperature and decreasing rainfall from northwest to southeast. Topography is rugged, especially in the northern half of the country, and most wild-land vegetation is evergreen, drought resistant and highly flammable. Forests and woodlands are a key element in the Portuguese landscape, covering more than one third of the country ([24]).

Fig. 1 - National Forest Inventory (NFI) plots used in the study. NFI 1995-1998 (n=102): black circles; NFI 2005-2006 (n=318): grey circles.

Results

Model fitting and selection

The wide range in biomass values observed for each shrub age reflects the high variability of the dataset used for modeling purposes, which encompass 420 forest plots (Fig. 2). Eucalyptus understory is dominated by Ulex spp., Erica spp. and Cistus ladanifer. In softwood stands the understory is characterized by the presence of Erica spp., Pterospartum tridentatum and Cistus ladanifer. In oaks plots, species of the genus Citysus are the most abundant. As for regeneration strategies, resprouter species were not present in 41% of the observations. A box-plot analysis revealed a lack of symmetry around the median (Fig. 2).

Fig. 2 - Comparative box-plots of shrub biomass (Mg ha-1) as a function of age (n=420). Grey circles and horizontal lines inside the box represent the mean and median, respectively. The upper and lower box limits represent the distance between the 25th and 75th percentiles and the whiskers extend to 10th and 90th percentile.

The general equation form selected for modeling the shrub biomass accumulation was as follows (eqn. 9):

(9)Biomi=(a+aresp+aG)⋅(−e−[bP+bresp+bslope+bT]t)" role="presentation"> Biomi=(a+aresp+aG)⋅(−e−[bP+bresp+bslope+bT]t)(9) (9)Biomi=(a1+a2resp+a3G)⋅(1−e−[b1P+b2resp+b3slope+b4T]t)

where Biom_i, is the shrub biomass (Mg ha^-1), resp is the resprouting percentage (R, %), G is stand basal area (m² ha^-1), P is precipitation (rain days year^-1), slope is in %, T is the mean annual temperature (°C), t is shrub age (years), and a_i and b_i are regression coefficients.

The five best-fitted equations are displayed in Tab. 5. All regression coefficients of the shrub biomass equations (a₁, a₂, a₃, b₁, b₂, b₃, b₄) were significant (p < 0.05) and biologically meaningful. Model Biom₁ had the best fit (model efficiency) and the smaller MSE, whereas models Biom₄ and Biom₅ showed the lowest R²_adj as well as the highest residual values (Tab. 5).

Tab. 5 - Adjusted coefficient of determination (R² _adj), mean square error (MSE) and parameter estimates (standard error in parentheses) of the regression coefficients (a_n, b_n) for the selected shrub biomass models (see eqn. 9). All the coefficients were significant at the 0.05 level.

Model	R 2 adj	MSE	Asymptote			Growth rate
			a 1	a 2	a 3	b 1	b 2	b 3	b 4
Biom 1	0.255	27.7	32.72 (6.47)	-0.239 (0.054)	-0.1528 (0.084)	-	0.00108 (0.00035)	0	0.00249 (0.00069)
Biom 2	0.252	28.74	36.55 (7.97)	-0.268 (0.066)	-0.1998 (0.0974)	0.00021 (0.00006)	0.00089 (0.00030)	0.00061 (0.00025)	0
Biom 3	0.227	29.72	31.77 (6.37)	-0.2286 (0.0639)	-0.159 (0.0897)	0.0003 (0.00008)	0.00109 (0.00036)	0	0
Biom 4	0.181	32.23	16.33 (2.41)	-0.0454 (0.0151)	-0.1582 (0.0666)	0.00107 (0.00023)	0	0	0
Biom 5	0.165	31.5	12.47 (1.42)	0	-0.115 (0.0635)	0.00133 (0.00025)	0	0	0

The maximum biomass accumulation, defined by the model asymptote (a_i parameters), is depending on stand variables, shrub regeneration strategies and the tree cover. No equation with a fixed asymptote was found for which model efficiency was acceptable. For all the five best-fitting equations, stand basal area (G, m² ha^-1) and/or the shrub biological traits affected the asymptote, i.e. the amount of shrub biomass decreased as G and the percentage of resprouters increased. Excluding G as a predictor and refitting the five models, the corresponding R²_adj and MSE for the different models were 0.24 and 31.0 (Biom₁), 0.23 and 32.0 (Biom₂), 0.21 and 32.4 (Biom₃), 0.17 and 35.4 (Biom₄) and 0.16 and 35.0 (Biom₅), respectively. Thus, the exclusion of G from the models decreased the prediction accuracy (R²_adj) and increased the errors (MSE).

Biomass growth rate (b_i parameters in eqn. 9) is affected by stand and shrub characteristics, as well as by topography and climate. Four out of the five selected equations included the precipitation, while the fifth included the temperature. As expected, higher precipitation and/or temperature increase the shrub growth rate, hence reaching the steady-state biomass faster. Prediction ability of the equations was considerably improved by including the resprouting-related variable. Only model Biom₂ included a topographic descriptor (the slope - Tab. 5). Elevation, aspect and the forest composition were not statistically significant in any of the model tested (p > 0.05).

Model validation and biological evaluation

All the five models described are ecologically meaningful, showing higher asymptotic values associated with lower canopy covers (as assessed by G) and quicker accumulation of shrubs when resprouters are present in the understory and for higher values of temperatures and/or precipitation. However, models Biom₄ and Biom₅ have a rather low asymptote and growth rate (Fig. 3). It is reasonable to expect slightly lower maximum biomass accumulation of shrubs in forest than in shrubland, but the two equations would considerably underestimate total biomass. They also present the lowest model efficiency and the least accurate estimates (i.e., high MSE). On the other hand, models Biom₁, Biom₂ and Biom₃ showed higher accuracy (highest R²_adj, lowest MSE - Tab. 5) and provided estimates of the maximum biomass accumulation in agreement with those reported in the literature (Fig. 3). Moreover, equations Biom₁ and Biom₄ had the lowest bias as revealed by the PRESS residuals (Tab. 6).

Fig. 3 - Shrub biomass predicted by the selected models as a function of shrub age using the following parameter values: precipitation = 103 days, temperature = 14.6 °C, resprouting =39.3%, basal area=5.65 m2 ha-1 and slope=14.6%.

Tab. 6 - Prediction Sum of Squares (PRESS), mean PRESS (MPRESS) and mean absolute PRESS residuals (MAPRESS) for the shrub biomass equations considered for the evaluation and validation.

Model	PRESS	MPRESS	MAPRESS
Biom 1	306.6	0.73	4.38
Biom 2	343	0.82	4.49
Biom 3	353.9	0.84	4.52
Biom 4	335.8	0.8	4.67
Biom 5	343.7	0.82	4.71

Model Biom₁ (which included resprouting percentage, basal area and temperature as independent variables) proved to be the most accurate model, while Biom₄ (one of the less biased models) was the less precise in terms of predictions (Tab. 6). To further characterize the selected Biom₁ model, graphs of predicted shrub biomass vs. observed biomass values (Fig. 4a) and vs. studentized residuals (Fig. 4b) were also plotted.

Fig. 4 - Predicted shrub biomass accumulation vs. observed biomass (a) and studentized residuals (b) in Portuguese forests using the model Biom1.

Model verification

The equation obtained for model Biom₁ closely describes the observed shrub accumulation pattern as a function of age (years), basal area and resprouter percentage (Tab. 4, Fig. 5). Despite the substantial variation in fuel accumulation, such equation may accurately predict the fuel build-up in various forest types. However, it should be noticed that the model was developed from data up to 15 years and its use for extrapolation on longer periods could be inappropriate, as suggested by the poor performances of the model at older ages. For extrapolation over periods longer than 15 years, the use of one of the difference equation forms seems more suitable (see below).

Fig. 5 - Comparative patterns of shrub biomass accumulation with stand age (after 15 years) and prescribed fire in a northeastern Portugal Pinus pinaster stands. Black symbols depict the measured shrub biomass and grey symbols represent the corresponding shrub biomass estimated by the model Biom1. Actual shrub biomass (Bio1 = 6.78 Mg ha-1, t1 = 10 years and Bio1 =10.48 Mg ha-1, t1 = 15 years) are connected by a dotted line to the value estimated by the difference eqn. 10 (Bio2). The solid biomass lines display the predicted (model Biom1) shrub biomass in Portuguese forests, the fixed values are the whole dataset means (resprouting = 15.4%, temperature= 11 °C).

Difference equation forms

Equation Biom₁ was used to obtain the difference equation form when age at time t₁ is known (eqn. 10):

(10)Bio=Bio⋅(a−aresp−aG)(−e−(bresp+bT)t)(a−aresp−aG)(−e−(bresp+bT)t)" role="presentation"> Bio=Bio⋅(a−aresp−aG)(−e−(bresp+bT)t)(a−aresp−aG)(−e−(bresp+bT)t)(10) (10)Bio2=Bio1⋅(a1−a2resp−a3G2)(1−e−(b2resp+b4T)t2)(a1−a2resp−a3G1)(1−e−(b2resp+b4T)t1)

where Bio₁ and Bio₂ are the biomass values (Mg ha^-1) at times t₁ and t₂, respectively; a₁ = 32.72; a₂ = 0.239; a₃ = 0.1528; b₂ = 0.00108; b₄ = 0.00249; G₁ and G₂ are the stand basal area (m² ha^-1) at t₁ and t₂, respectively; resp is the percentage (%) of resprouters, assumed as fixed over time (15.4%), and T is the mean annual temperature (11 °C).

In addition, rearranging eqn. 9 the above model can also be used when age (t₁) is not known, making it more useful for practical applications. Solving the equation Biom₁ for t, and substituting this expression into eqn. 11, the correspondent difference equation form is derived as (eqn. 11):

(11)Biot+i=a[−e−(bresp+bT)−(bresp+bT)ln⁡(−Biota)+i]" role="presentation"> Biot+i=a[−e−(bresp+bT)−(bresp+bT)ln(−Biota)+i](11) (11)Biot+i=a[1−e−(b2resp+b4T)1−(b2resp+b4T)ln⁡(1−Biota)+i]

where Bio_t and Bio_t+í are the biomass (Mg ha^-1) at ages t and t+i, respectively; a = 32.72 - 0.239 resp - 0.1528 G; b₂ = 0.00108; b₄ = 0.00249; i is the projection length; G is the stand basal area (m² ha^-1) at t; resp is the percentage (%) of resprouters assumed as invariable over time, and T is the mean annual temperature.

Biomass accumulation may vary for several reasons (e.g., differences in site productivity, fire severity, climate, and competition) and part of the relevant factors may be unknown. The difference equation (eqn. 10) allows the estimation of future biomass at time t₂ when biomass at time t₁ is known (i.e., information on the age and the shrub biomass at t₁). For demonstration purposes, based on the information displayed in Fig. 5, the actual shrub load from the maritime pine stand chronosequence was used to estimate the biomass Bio₂ at time t₂, i.e., the initial shrub load (Bio₁ =6.74 Mg ha^-1 and 10.48 Mg ha^-1) at time t₁ (10 and 15 years) was taken as the basis to predict biomass Bio₂ at time t₂ (13 and 18 years), with a stand basal area G₁ (32.1 and 6.9 m² ha^-1) and G₂ (19.7 and 15 m² ha^-1), respectively. For more details, the values of parameters describing the Fig. 5 are listed in Tab. 4. The applicability of the eqn. 10 is clearly depicted in Fig. 5 (black dotted lines for the range t₁ = 10 years to t₂ = 13 years and grey dotted lines for the range t₁ = 15 years to t₂ = 18 years). Indeed, the estimated values Bio₂ (8.86 Mg ha^-1 and 11.34 Mg ha^-1, eqn. 10) are similar to the measured biomass (black symbol in Fig. 5). For shrub ages below 15 years the predictions of Biom₁ model (grey symbol in Fig. 5) closely match the observed shrub accumulation. Therefore, eqn. 10 predicts quite well future biomass for accumulation periods above 15 years.

If information about initial shrub age (t₁) is missing, the difference equation form (eqn. 11) should be used to project the future shrub biomass, for it does not use age as an explicit variable. In these cases we do not assume that the shrubs follow strictly the curves of biomass accumulation prediction (model Biom₁) as it may start from a different level of biomass.

Thereby, the use of the difference equations forms implies one of the following situations: (i) when only one measurement is available and age is known at t₁, biomass can be obtained from eqn. 10; (ii) when only one measurement is available and age is unknown, with no other additional information, the biomass may be estimated by eqn. 11. If no information on shrub load for a specific age is available, the model Biom₁ should be used for initialization and prediction of forest shrub accumulation.

Implications for forest management - application examples

For forest management purposes, the model Biom₁ was used to compute the shrub biomass in a hypothetical stand using the mean values of the independent variables considered (Fig. 6). Basal area, used as an indicator of stand-level competition, ranged from 5.5 m² ha^-1 to 49.5 m² ha^-1. In general, an increase in the basal area implies a higher canopy cover, reducing light and water availability for the shrub layer. Thus, increasing competition for resources reduces shrub growth over time (Fig. 6). Setting t = 15, a difference in fuel load of about 5 Mg ha^-1 was observed for extreme basal area values.

Fig. 6 - Shrub biomass predicted over time as a function of the basal area using the model Biom1 for a stand characterized by mean conditions (temperature = 14.6 °C, resprouting = 39.3%).

Fig. 7 displays how shrub biomass is impacted by shrub regeneration strategy (two levels of resprouter abundance) and by basal area (two extreme values of basal area) using model Biom₁. The lowest biomass usually appears in the resprouters’ community (R = 100%), particularly in the highest basal area stands (Fig. 7a). In contrast, biomass recovery is faster and reaches higher levels when the understory is comprised exclusively of seeders (i.e., R = 0%), particularly at the lowest basal area level, where 13-15 Mg ha^-1 is reached after 15 years (Fig. 7b).

Fig. 7 - Shrub biomass over time predicted by the model Biom1 using two levels of resprouter abundance (0 and 100%) and two levels of stand basal area (49.5 and 5.5 m2 ha-1).

Discussion and Conclusions

Mediterranean-type shrublands and forests are highly variable in aboveground biomass ([40]). Biomass estimation in different forest stands is a difficult task due to structural heterogeneity and the dynamic nature of vegetation. Given current knowledge gaps in this topic, the main objective of the present study was to model shrub accumulation over time in the presence of a tree overstorey. Our starting hypothesis was that shrub growth under a tree canopy should depend on both stand variables and climate, reflecting their importance to biomass growth. Prediction models were developed based on NFI data of diverse geographical origin and various floristic compositions, in order to extend their applicability to most Portuguese forests with an abundant understory vegetation, especially those dominated by Erica spp., Ulex spp., Cytisus spp. or Cistus ladanifer. The strategy adopted can be viewed as an advantage over site-specific models ([18]).

Dynamics of shrub biomass under tree cover in the Mediterranean climate, and its relation to stand variables had been poorly studied. This study extends previous investigations by introducing measurable stand variables and shrub regenerating strategy in a time-dependent shrub biomass model, using a non-linear regression technique based on Olson’s model. Moreover, two alternative difference equation forms obtained directly from the selected model could be used to project the observed shrub biomass: difference eqn. 10 when age is known and difference eqn. 11 when age is not available. This will increase the quality of the shrub build-up models when shrub age is not readily available. Moreover, to provide information on biomass trends over time, the proportion of resprouters in the understorey is the sole shrub information required by our model. Indeed, field assessment of shrub regeneration strategy has the advantage of being relatively straightforward.

Among the studied models, Biom₁ provided the best performances for shrub biomass estimation under forest canopies in Portugal. Also models Biom₄ and Biom₅ provided biologically meaningful predictions, but their asymptote and growth rate parameters are low and inconsistent with the range of values reported in the literature ([65]). On one hand, the function Biom₁ depends on a shrub trait (regeneration mode), affecting both growth rate and maximum biomass, and on a descriptor of the tree overstorey (stand basal area), an indicator of stand-level competition. On the other hand, climatic conditions such as temperature influence vegetative growth rate. As expected, forest shrub biomass estimated with our models is lower than that reported for shrublands in the Mediterranean (e.g., [75]). Furthermore, our estimates of shrub biomass grown under tree canopy are in line with values reported by previous studies in Portugal, which do not exceed 24 Mg ha^-1 ([27]).

Rosa et al. ([67]) indicates that shrublands in Portugal require 15 years to reach a steady-state aboveground biomass. Our results indicate that longer periods are needed when tree cover is present. Obviously, such slower build-up of shrub biomass may be explained by the competition for available resources between the understory and the overstory ([43]). Similarly, Castedo-Dorado et al. ([14]) showed that selected overstory variables could be suitable indicators of the relative availability of light, nutrients, water or growing space in modeling the maximum shrub development. Nevertheless, several authors reported that the amount of understory vegetation in maritime pine and eucalyptus stands is practically the same 10 and 30 years after fire. Indeed, after an initial period of relatively vigorous growth, understory biomass tends to stabilize when shrubs reach their adult size ([73]).

Shrub biomass tends to increase with time since disturbance, with variations depending on the shrub type. However, the highest accumulation usually occurs in seeders communities, particularly when basal area is low (Fig. 7). This agrees with the general knowledge on different growth patterns between resprouters and non-resprouters. Non-resprouters begin flowering earlier and more abundantly after disturbance ([58], [5]), allowing a faster colonization of suitable microhabitats as compared to resprouting species, whose regrowth from established root systems is restricted to microsites previously occupied ([41], [11]). Moreover, resprouters need more energy and time to regenerate after disturbance, due to their resource allocation to the replacement of damage tissues ([5], [2]). According to our results, several authors (e.g. [5], [51]) reported that many resprouting species have lower growth rates than non-resprouters of comparable age (Fig. 7). In Mediterranean-type ecosystems, frequent fires decrease the abundance of seeders, while at intermediate fire-return intervals, they are favored in comparison to resprouters ([76], [59], [45]).

Overstory basal area in the proposed models is a proxy for competition, decreasing shrub biomass with higher tree stocking. This general trend is consistent with other studies, reporting that the maximum shrub development was also limited by overstory variables such as basal area (G - [21], [14]). Thus G serves as an indicator of stand competition, and has the advantage of being relatively simple to obtain in the field or to be inferred from growth and yield models ([14]).

In a Mediterranean-type climate, the regional patterns of vegetation structure and composition are mainly dependent on temperature and water availability ([44], [40]). Hence it was not surprising to find temperature as a significant predictor in proposed model. Most phenological models include temperature as a proxy for developmental rates ([16]). For instance, higher air temperature values in spring strongly induces an earlier start of plant development ([22], [46], [17]). Moreover, temperature directly affects soil moisture, which is the main source of water for plant growth ([44]).

Forest composition did not contribute to explain variation in shrub biomass. However, post fire composition, structure and richness of the whole plant community are directly related to the re-establishment of the dominant species in the canopy ([45], [9]). Differences in shrub biomass accumulation among the different forest types considered in this study were generally poor; Maritime pine plots attained higher values and maximum biomass was qualitatively similar to other species. It is important to remark that forest composition in Portugal is strongly determined by human action and it is often decoupled from the potential vegetation type. The dominant shrub communities generally occur independently of overstorey composition. Consequently, variation in basal area and site quality are likely to overwhelm the influence of forest composition on shrub biomass.

Several limitations of the proposed models have also to be highlighted here. The variability of the initial dataset was considerably high, lowering the model efficiency. Variation in shrub species composition, the use of generic bulk densities to estimate biomass, and site-specific factors are certainly involved in the modeling performance. Similarly, several models with low R² values (indicating that some variation remains unexplained) were considered useful to support optimal forest management decisions, for instance, mushrooms yields in pine forest planning ([6], [56]). In any case, we consider this study to be an important step forward, since similar models have seldom been developed for Europe.

Future research would greatly benefit from the existence of information from permanent shrub/biomass plots measured over time, in order to better understand biomass dynamics (e.g., shrub senescence) and carbon sequestration rates at different spatial scales. Additionally, variables such as soil type and summer water deficit may contribute to improve the model performances.

The successful management of fuel load in fire-prone regions is a challenging task that calls for the integration of forest and fire management activities in order to decrease fire hazard. Several studies indicate that fuel treatments (i.e., reduction of fuels in forests) change wildfire behavior and severity and enhance the effectiveness of fire suppression operations ([50]). Models such those developed in this study have practical applications in the assessment of fire hazard and in the definition of general prescriptions to plan fuel treatments. In fact, such models allow to quantify the impact of silviculture operations and help to define management options that may decrease wildfire occurrence. For instance, a preliminary shrub build-up model ([7]) was integrated with a growth and yield model to accomplish maritime pine stand-level optimizations and determine optimal stand-level treatments (e.g., thinning, fuel treatment), so as to reduce the hazard of fire ([32]). This information is very valuable as it may effectively support the development of adaptive management strategies ([31], [32]).

Investigations on carbon flux implications of fuels reduction treatments are of increasing interest ([39]). Fire-related climate change mitigation options include decreasing emissions through fuel reduction treatments and using the removed biomass to replace fossil fuels for energy production ([13], [47]). The shrub build-up model can be used by forest managers to predict fuel loads in the frame of understory removal to decrease fire hazard. It may also have practical application for ecologists, allowing the estimation of carbon storage in the understory and the assessment of how future wildfire emissions will change in response to fuel treatments, helping to reduce the uncertainty in emission estimates. This analysis may be done in two steps. If no information on shrub biomass is available the proposed model would be used for initialization (i.e. estimate initial shrub biomass) and to predict curves of forest shrub accumulation without previous information on shrub load over a specific time span; otherwise one of the difference equation forms are recommended to estimate future biomass from the initial biomass.

Despite disregarding site-specific conditions in relation to shrub composition and some simplifications in regards to the whole fuel complex, it is reasonable to conclude that the combination of the proposed model and the difference equation forms is able to adequately describe the typical shrub biomass accumulation in Portuguese forests (Fig. 5). Additionally, our biomass models provide sound estimates of biomass growth on the short term. For longer periods (more than 15 years), information on regeneration, mortality, thinning and succession has to be taken into consideration in order to increase the accuracy of biomass estimates. These general equations are expected to help forest management decision-making as a tool to support decisions on where and when fuel treatments are required.

Acknowledgments

This research was conducted in the frame of the project PTDC/AGR-CFL/64146/2006 “Decision Support Tools for Integrating Fire and Forest Management Planning” and project FIRE-ENGINE “Flexible Design of Forest Fire Management Systems” (MIT/ FSE/ 0064/2009), funded by the Portuguese Science Foundation (FCT), and and partly funded by the project ForEAdapt “Knowledge exchange between Europe and America on forest growth models and optimization for adaptive forestry”, under grant agreement no. PIRSES-GA-2010-269257 and INTEGRAL “Future Oriented Integrated Management of European Forest Lands, both funded by the European Union Seventh Framework Programme (FP7-PEOPLE-2010 -IRSES). The FCT support for funding the doctoral program plan of Brigite Roxo Botequim (SFRH-BD-44830-2008) and Susete Marques (SFRH/BD/62847/2009) is also acknowledged. JGG participated in this research under the framework of the Project PTDC/AGR-FOR/4526/2012 “Models and Decision Support Systems for Addressing Risk and Uncertainty in Forest Planning” (SADRI).

Authors would also like to thank the financial support for a post-fellowship by the FCT (SFRH/BPD/63979/ 2009) and by the University of Eastern Finland. Finally, the Portuguese Forest Service (ICNF) is acknowledged for supplying the NFI databases.

References

(1) Amaro A, Reed D, Tomé M, Themido I (1998). Modeling dominant height growth: eucalyptus plantation in Portugal. Forest Science 44: 37-46. Online | Gscholar

(2) Arnan X, Rodrigo A, Retana J (2007). Post-fire regeneration of Mediterranean plant communities at a regional scales is dependent on vegetation type and dryness. Journal of Vegetation species 18: 11-122.

(3) Baeza M, Raventos J, Escarré A, Vallejo VR (2006). Fire risk and vegetation structural dynamics in Mediterranean shrubland. Plant Ecology 187 (2): 189-201.

(4) Brown JK (1971). A planar intersect method for sampling fuel volume and surface area. Forest Science 17 (1): 96-102. Online | Gscholar

(5) Bond WJ, Van Wilgen BW (1996). Fire and plants (1st edn). Population and Community Biology Series 14, Chapman and Hall, London, UK, pp. 263.

(6) Bonet JA, Pukkala T, Fischer CR, Palahí M, Martinez de Aragón J, Colinas C (2008). Empirical models for predicting the production of wild mushrooms in Scots pine (Pinus sylvestris L.) forests in the Central Pyrenées. Annals of Forest Science 65 (2): 206-206.

(7) Botequim B, Borges P, Carreiras J, Oliveira MM, Borges J (2009). Development of a shrub growth model in understory conditions (preliminary model). Technical Report 7, FORCHANGE, Instituto Superior de Agronomia, Lisboa, Portugal, pp. 12. Gscholar

(8) Botequim B, Garcia-Gonzalo J, Marques S, Ricardo A, Borges JG, Oliveira MM, Tomé J, Tomé M (2013). Assessing wildfire risk probability in Eucalyptus globulus Labill stands in Portugal. iForest 6: 217-227.

(9) Broncano MJ, Retana J, Rodrigo A (2005). Predicting the recovery of Pinus halepensis and Quercus ilex forests after a large wildfire in northeastern Spain. Plant Ecology 180: 47-56.

(10) Burkhart HE, Tomé M (2012). Modeling forest trees and stands. Springer, Berlin, Germany, pp. 457.

(11) Calvo L, Tárrega R, Luis E, Valbuena L, Marcos E (2005). Recovery after experimental cutting and burning in three shrub communities with different dominant species. Plant Ecology 175-185.

(12) Canadell J, Zedler PH (1995). Underground structures of woody plants in Mediterranean ecosystems of Australia, California, and Chile. In: “Ecology and Biogeography of Mediterranean Ecosystems in Chile, California and Australia” (Arroyo MTK, Zedler PH, Fox Deds M). Springer-Verlag, New York, USA, pp. 177-210.

(13) Canadell JG, Raupach MR (2008). Managing forests for climate change mitigation. Science 320: 1456-1457.

(14) Castedo-Dorado F, Gómez-Vásquez I, Fernandes PM, Crecente-Campo F (2012). Shrub fuel characteristics estimated from overstory variables in NW Spain pine stands. Forest Ecology and Management 275: 130-141.

(15) Castro H, Freitas H (2009). Above-ground biomass and productivity in the Montado: from herbaceous to shrub dominated communities. Journal of Arid Environments 73 (4-5): 506-511.

(16) Cesaraccio C, Spano D, Duce P, Snyder RL (2001). An improved model for determining degree-day values from daily temperature data. International Journal of Biometeorology 45 (4): 161-169.

(17) Chmielewski FM, Rötzer T (2001). Response of tree phenology to climate change across Europe. Agricultural and Forest Meteorology 108 (2): 101-112.

(18) Chojnacky DC (1992). Estimating volume and biomass for dryland oak species. In: Proceedings of the Conference “Ecology and Management of oaks and Associated Woodlands: Perspectives in the Southwestern United States and Northern Mexico”. Sierra Vista (Arizona, USA) 27-30 April 1992. Rocky Mountain Forest and Range Experiment Station, USDA Forest Service, Fort Collins, CO, USA, pp. 151-161. Online | Gscholar

(19) Cieszewski CJ, Bailey RL (2000). Generalized algebraic difference approach: theory based derivation of dynamic equations with polymorphism and variable asymptotes. Forest Science 46: 116-126. Gscholar

(20) Clutter JL, Fortson JC, Pienaar LV, Brister GH, Bailey RL (1983). Timber management: a quantitative approach. John Wiley and Sons, Inc, New York, USA, pp. 333. Gscholar

(21) Coll L, González-Olabarria JR, Mola-Yudego B, Pukkala T, Messier C (2011). Predicting understory maximum shrubs cover using altitude and overstory basal area indifferent Mediterranean forests. European Journal of Forest Research 130: 55-65

(22) Correia OA, Martins AC, Catarino FM (1992). Comparative phenology and seasonal foliar nitrogen variation in Mediterranean species of Portugal. Ecologia Mediterranea 18: 7-18. Online | Gscholar

(23) DGF (2001). Inventário Florestal Nacional: Portugal Continental (3ª Revisão) [Forest National Inventory: Portugal (3rd revision)]. Direcção Geral das Florestas, Lisboa, Portugal, pp. 233. [in Portuguese] Gscholar

(24) DGRF (2006). Resultados do Inventário Florestal Nacional 2005/2006 (5ª Revisão) [Results from the National Forest Inventory 2005/2006 (5th Revision)]. Direcção Geral dos Recursos Florestais, Lisboa, Portugal, pp. 70. [in Portuguese] Gscholar

(25) Fernandes P, Rego FC (1998). Equations for estimating fuel load in shrub communities dominated by Chamaespartium tridentatum and Erica umbellata. In: Proceedings of the “3rd International Conference on Forest Fire Research & 14th Fire and Forest Meteorology Conference” (Viegas DX ed). Luso (Portugal) 16-20 Nov 1998. ADAI, University of Coimbra, Coimbra, Portugal, pp. 2553-2564. Gscholar

(26) Fernandes PA, Loureiro CA, Botelho HS (2004). Fire behaviour and severity in a maritime pine stand under differing fuel conditions. Annals of Forest Science 61 (6): 537-544.

(27) Fernandes PM, Rigolot E (2007). The fire ecology and management of maritime pine (Pinus pinaster Ait.). Forest Ecology and Management 241 (1-3): 1-13.

(28) Fernandes PM (2009a). Combining forest structure data and fuel modelling to classify fire hazard in Portugal. Annals of Forest Science 66 (4): 415-415.

(29) Fernandes PM (2009b). Examining fuel treatment longevity through experimental and simulated surface fire behaviour: a maritime pine case study. Canadian Journal of Forest Research 39 (12): 2529-2535.

(30) Fernández Alés R, Martín A, Ortega F, Alés EE (1992). Recent changes in landscape structure and function in a Mediterranean region of Spain (1950-1984). Landscape Ecology 7: 3-18.

(31) Ferreira L, Constantino M, Borges J, Garcia-Gonzalo J (2012). A stochastic dynamic programming approach to optimize short-rotation coppice systems management scheduling: an application to eucalypt plantations under wildfire risk in Portugal. Forest Science 58 (4): 353-365.

(32) Ferreira L, Constantino M, Borges JG (2014). A stochastic approach to optimize Maritime pine (Pinus pinaster Ait.) stand management scheduling under fire risk: an application in Portugal. Annals of Operations Research 219 (1): 359-377.

(33) Garcia-Gonzalo J, Marques S, Borges JG, Botequim B, Oliveira MM, Tomé J, Tomé M (2011). A three-step approach to post-fire mortality modelling in Maritime pine (Pinus pinaster Ait) stands for enhanced forest planning in Portugal. Forestry 84 (2): 197-206.1

(34) Garcia-Gonzalo J, Zubizarreta-Gerendiain A, Ricardo A, Marques S, Botequim B, Borges J. G, Oliveira MM, Tomé M, Pereira JMC (2012). Modelling wildfire risk in pure and mixed forest stands in Portugal. Allgemeine Forst und Jagdzeitung (AFJZ) - German Journal of Forest Research 183 (11/12): 238-248. Gscholar

(35) Garcia-Gonzalo J, Pukkala T, Borges J (2014). Integrating fire risk in stand management scheduling. An application to Maritime pine stands in Portugal. Annals of Operations Research 219 (1): 379-395.

(36) González JR, Palahi M, Pukkala T (2006). A fire probability model for forest stands in Catalonia (north-east Spain). Annals of Forest Science 63: 1-8.

(37) Gould JS, McCaw WL, Cheney NP (2011). Quantifying fine fuel dynamics and structure in dry eucalypt forest (Eucalyptus marginata) in Western Australia for fire management. Forest Ecology and Management 262 (3): 531-546.

(38) Hosmer DW, Lemeshow S (2000). Applied logistic regression (2nd edn). Series “Probability and Mathematical Statistics”, no. 452, Wiley, New York, USA, pp. 307. Gscholar

(39) Hurteau MD, North M (2010). Carbon recovery rates following different wildfire risk mitigation treatments. Forest Ecology and Management 260: 930-937.

(40) Joffre R, Rambal S (2002). Mediterranean Ecosystems. In: “Encyclopedia of Life Science” (edn M.P. Ltd). Nature Publishing Group, London, UK, pp. 1-7. Gscholar

(41) Keeley JE, Zedler PH (1978). Reproduction of Chaparral shrubs after fire: a comparison of sprouting and seeding strategies. American Midland Naturalist 99: 142-161.

(42) Keeley JE (1986). Resilience of Mediterranean shrub communities to fire. In: “Resilience in Mediterranean-type ecosystems” (Dell B, Hopkins AJM, Lamont Beds B). Dr W Junk Publishers, Dordrecht, The Netherlands, pp. 95-112. Gscholar

(43) Kozlowski TT, Kramer PJ, Pallardy SG (1991). The physiological ecology of woody plants. Academic Press Inc, San Diego, CA, USA, pp. 30. Gscholar

(44) Kummerow J, Montenegro G, Krause D (1981). Biomass, phenology and growth. In: “Resource Use by Chaparral and Matorral” (Miller PC ed). Ecological studies, vol. 39, Springer-Verlag, New York, USA, pp. 69-96.

(45) Lloret F, Pausas J G, Vila M (2003). Responses of Mediterranean plant species to different fire frequencies in Garraf Natural Park (Catalonia, Spain): field observations and modelling predictions. Plant Ecology 167: 223-235.

(46) Maak K, von Storch H (1997). Statistical downscaling of monthly mean air temperature to the beginning of flowering of Galanthus nivalis L. in northern Germany. International Journal of Biometeorology 41 (1): 5-12.

(47) Malmsheimer RW, Heffernan P, Brink S, Crandall D, Deneke F, Galik C, Gee E, Helms JA, McClure N, Mortimer M, Ruddell S, Smith M, Stewart J (2008). Forest management solutions for mitigating climate change in the United States. Journal of Forestry 106 (3): 115-171. Gscholar

(48) Marques S, Garcia-Gonzalo J, Botequim B, Ricardo A, Borges JG, Tomé M, Oliveira MM (2012). Assessing wildfire risk probability in Pinus pinaster Ait. stands in Portugal. Forest Systems 21 (1): 111-120.

(49) McCaw L, Neal J, Smith R (2002). Stand characteristics and fuel accumulation in a sequence of even-aged karri (Eucalyptus diversicolor) stands in south-west Western Australia. Forest Ecology and Management 158: 263-271.

(50) Mercer DE, Haigth RG, Prestemon JP (2008). Analyzing trade-offs between fuels management, suppression, and damages from wildfire. In: “The Economics of Forest Disturbances: Management, Suppression, and Damages from Wildfire” (Holmes TP, Prestemon JP, Abt KL eds). Series “Forestry Sciences”, vol. 79, Springer, The Netherlands, pp. 247-272.

(51) Midgley JJ (1996). Why the world’s vegetation is not totally dominated by resprouting plants; because reprouters are shorter than reseeders. Ecography 19: 92-95.

(52) Myers RH (1990). Classical and modern regression with applications. PWS-Kent Publishing, Boston, USA, pp. 488. Gscholar

(53) Navar J, Nájera J, Jurado E (2001). Preliminary estimates of biomass growth in the Tamaulipan thornscrub in north-eastern Mexico. Journal of Arid Environments 47: 281-290.

(54) Olson JS (1963). Energy storage and the balance of producers and decomposers in ecological systems. Ecology 44: 322-331.

(55) Palahí M, Tomé M, Pukkala T, Trasobares A (2004). Site index model for Pinus sylvestris in north-east Spain. Forest Ecology and Management 187: 34-47.

(56) Palahí M, Pukkala T, Bonet JA, Colinas C, Fischer CR, Martínez de Aragón JR (2009). Effect of the inclusion of mushroom values on the optimal management of even-aged Pine stands of Catalonia. Forest Science 55 (60): 503-511(9). Online | Gscholar

(57) Parresol BR (1999). Assessing tree and stand biomass: a review with examples and critical comparisons. Forest Science 45 (4): 573-593. Online | Gscholar

(58) Pate JS, Froend RH, Bowen BJ, Hansen A, Kuo J (1990). Seedling growth and storage characteristics of seeder and resprouter species of Mediterranean-Type ecosystems of SW Australia. Annals of Botany 65: 585-601. Online | Gscholar

(59) Pausas JG (2001). Resprouting vs. seeding - a Mediterranean perspective. Oikos 94: 193-194.

(60) Pausas JG (2004). Changes in fire and climate in the eastern Iberian Peninsula (Mediterranean basin). Climatic Change 63: 337-350.

(61) Pausas JG, Bradstock RA, Keith DA, Keeley JE (2004). The GCTE Network. Plant functional traits in relation to fire in crown-fire ecosystems. Ecology 85: 1085-1100.

(62) Pausas JG, Keeley JE, Verdú M (2006). Inferring differential evolutionary processes of plant persistence traits in Northern Hemisphere Mediterranean fire-prone ecosystems. Journal of Ecology 94: 31-39.

(63) Pereira JMC, Santos TN (2003). Fire risk and burned area mapping in Portugal. Direcção Geral das Florestas, Lisboa, Portugal. [in Portuguese] Gscholar

(64) Plucinski MP (2003). The investigation of factors governing ignition and development of fires in heathland vegetation. PhD thesis, University of New South Wales - Australian Defence Force Academy, School of Mathematics and Statistics, Sidney, Australia, pp. 366. Gscholar

(65) Rambal S, Hoff C (2001). Mediterranean ecosystems and fire: the threats of global change. In: “Large Forest Fires” (Moreno JM ed). Backbuys Publishers, Leiden, The Netherlands, pp. 187-213. Gscholar

(66) Richards FJ (1959). A flexible growth function for empirical use. Journal of Experimental Botany 10 (2): 290-301.

(67) Rosa IM, Pereira JM, Tarantola S (2011). Atmospheric emissions from vegetation fires in Portugal (1990-2008): estimates, uncertainty analysis, and sensitivity analysis. Atmospheric Chemistry and Physics 11: 2625-2640.

(68) SAS Institute Inc (2000). SAS/STAT user’s guide (v. 8 edn). SAS Institute Inc, Cary, NC, USA. Gscholar

(69) Schmidt KM, Menakis JP, Hardy CC, Hann WJ, Bunnell DL (2002). Development of coarse-scale spatial data for wildland fire and fuel management. Gen. Tech. Rep. RMRS-GTR-87, Rocky Mountain Research Station, USDA Forest Service, Fort Collins, CO, USA, pp. 41. Gscholar

(70) Schumacher FX (1939). A new growth curve and its application to timber-yield studies. Journal of Forestry 37: 819-820. Gscholar

(71) Tomé M, Oliveira T, Soares P (2006a). O modelo Globulus 3.0. [The Globulus model 3.0.] Publicações GIMREF-RC2/2006, Universidade Técnica de Lisboa. Instituto Superior de Agronomia. Centro de Estudos Florestais. Lisboa, Portugal, pp. 23. Gscholar

(72) Tomé J, Tomé M, Barreiro S, Paulo JA (2006b). Age-independent difference equations for modelling tree and growth. Canadian Journal of Forest Research 36: 1621-1630.

(73) Trabaud L, Grosman J, Walter T (1985). Recovery of burnt Pinus halepensis Mill. Forest. I. Understorey and litter phytomass development after wildfire. Forest Ecology and Management 12 (3-4): 269-277.

(74) Vanclay JK, Skovsgaard JP (1997). Evaluating forest growth models. Ecological Modelling 98:1-12.

(75) Vega JA, Fernandes P, Cuinas P, Fonturbel M, Perez J, Loureiro C (2006). Fire spread analysis of early summer field experiments in shrubland fuel types of northwestern Iberia. Forest Ecology and Management 234: S102. Gscholar

(76) Zedler PH (1995). Fire frequency in southern California shrublands: biological effects and management options. In: “Brushfires in California wildlands: ecology and resource management” (Keeley JE, Scott T). International Association of Wildland Fire, Fairfield, Washington, DC, USA, pp. 101-112. Gscholar