Chapter 18 Forest dynamics model
This chapter provides an overview of a forest dynamics model, which builds on the previous models and allows simulating the recruitment, growth, mortality and management of woody plant cohorts in a forest stand. The model was described in De Cáceres et al. (2023) and is run using function fordyn() for a set of years.
18.1 Design principles
The design of the forest dynamics model is to a large degree inherited from the water balance and growth models it builds on. Readers should refer to former sections to learn the design of the basic water balance (see 3.1), advanced water/energy balance (see 7.2) or growth/mortality processes (see 15.1). Regeneration, from either seed germination and establishment or via resprouting, and forest management are explicitly simulated at the level of fordyn() and hence will be the focus of this chapter.
18.1.1 Recruitment from seeds
Recruitment of saplings involves a number of processes (flowering and pollination, fruit/seed production, dispersal, storage, seed predation, germination, seedling establishment and survival until the sapling stage). All these processes have their own biotic and abiotic drivers (Price et al. 2001), so that modelling becomes a challenging task. Processes leading to sapling recruitment are frequently extremely simplified or their mechanisms ignored in many forest models (Price et al. 2001). The design of the forest dynamics model with respect to recruitment follows that of many gap models. Local seed production is considered in a binary way, where plants are considered fertile and able to produce viable seeds if they reach a given height (different for shrubs and trees). Alternatively, the user can specify a set of species whose seeds arrive to the target stand via simulation control parameters or dispersal processes (see chapter 20). A seed bank exists in the forest stand, where seed levels are relative (up to 100%) and result from the interplay between seed rain and seed mortality. Actual recruitment from seeds depends on a set of regeneration thresholds are used to determine whether recruitment of new saplings (i.e. ingrowth into a given diameter class for trees or size for shrubs) is possible. Typically, regeneration thresholds concern environmental conditions, although some models also consider ungulate browsing (Wehrli et al. 2007). In our case we focus on three environmental drivers limiting the transition from seedlings to saplings:
- Tolerance to low temperatures, indicated by the mean temperature of the coldest month.
- Drought tolerance, indicated by the annual moisture index (annual precipitation divided by annual evapotranspiration).
- Shade tolerance, indicated by the percentage of photosynthetic active radiation reaching the ground.
A target species can recruit on a given year if current conditions for all three environmental indicators are above the species tolerance thresholds. A constant probability of recruitment determines actual recruitment within these bioclimatic limits. Recruitment densities and plant size of recruited individuals are specified via species parameters. Trees recruited from seed are subject to self-thinning processes before attaining the diameter of ingrowth (see 17.4.1).
18.1.2 Resprouting
Resprouting is a common feature of Mediterranean species. Most, if not all, Mediterranean broadleaved species are able to resprout from buds protected in branches, subterranean structures (e.g. burl, taproot, lateral roots) or in the root collar after a disturbance destroys the aerial part. Yet, some differences in the resprouting ability (i.e. survival) exist regarding the type of disturbance and the species affected.
In medfate, resprouting is assumed happen from buds at the stump/stool level after a given disturbance has led to aboveground plant die-back. Resprouting can occur as a response to different disturbances. For example, resprouting occurs after clipping (e.g. if management is applied), after a fire, or after cavitation has led to dessication of above-ground plant organs.
Regarding disturbance type, differences in survivorship through resprouting mostly depends on the intensity of the disturbance experienced, with an increasing occurrence of mortality after disturbances in the following order: browsing, cutting, drought and fire (Espelta et al. 1999). Fire is the disturbance that causes greatest mortality, probably because it may physically destroy part of the bud bank due to lethal temperatures or directly charring it (Pausas et al. 2016). Moreover, resprouting vigor also depends on the size of the plant (Moreira et al. 2012). According to different literature sources, the number of resprouts initially produced is approximately 2 per 1 cm2 of stump area (Retana et al. 1992; García‐Jiménez et al. 2017). Hence, in medfate the final density of resprouts depends on the disturbance that caused resprouting (via survivorship) and on the diameter of the parent plant, but in all cases the resprouts inherit the root system of their parent.
Like saplings originated from seeds, resprouts are subject to self-thinning processes before attaining the diameter of ingrowth (see 17.4.1).
18.1.3 Forest management
Forest management is an optional process in fordyn() simulations. Furthermore, management actions need to be defined in an external function to be supplied by the user. This design was chosen because there are multiple potential management strategies, and different traditions are followed in different countries. To work properly within fordyn() simulations, the supplied management function needs to accept a forest object as input, as well as a list of management parameters, and it has to return the reduction in tree density, the reduction in shrub cover and the management parameters to be applied in subsequent calls.
Although the former design allows users to tailor management functions to their simulation needs, the package includes defaultManagementFunction() to facilitate the simulation of these process in many situations, and which is briefly described here. The function implements two different management models:
- Irregular management: An uneven-aged stand is managed using thinning operations each time a threshold value of a chosen stand-level metric (such as basal area or Hart-Becking index) is trespassed. Thinning operations can focus on trees of specific diameter classes.
- Regular management: An even-aged (monospecific) stand is managed in cycles where thinning (preparatory) cuts are followed by final (regeneration) cuts. Thinning operations are similar as those of irregular management. One or several final cuts can be scheduled, the first starting whenever mean diameter surpasses a chosen threshold, and the following after a chosen number of years. Optionally, tree planting (of a chosen species) can be scheduled to occur after the last final cut.
The former models apply cuts on tree cohorts irrespective of their species. The default management function applies shrub clearing each time there is an operation on trees, resulting in removal of any shrub cover above a given maximum value.
18.2 State variables
The main state variables of the forest dynamic model are those conforming the structure and composition of the forest stand, i.e. the set of woody cohorts (either trees or shrubs) and their attributes (height, density, DBH, cover, etc.). Additionally, the seed species identity and relative abundance in the seed bank are also state variables of the model. Since the model performs calls to the growth() model, many other state variables are defined for intra-annual simulations (see 15.4). When including forest management action, additional state variables may be defined as management parameters.
18.3 Process scheduling
The fordyn() model divides the period to be simulated in years, which is the top-level time step of simulations. Given an input forest object, the function first initializes the input for function growth(). For each year to be simulated the model the performs the following steps:
- Calls function
growth()to simulate daily water/carbon balance, growth and mortality processes (sub-daily processes may also be involved intranspirationMode = "Sperry"ortranspirationMode = "Sureau"). See section 15.4 for details of growth scheduling. - If a management function is supplied, calls this function and apply the resulting reductions in tree density and shrub cover.
- If required, simulates seed production, seed bank mortality, seed rain, recruitment from seeds and/or resprouting. In a spatial context, seed rain will include seeds dispersed from other forest stands 20.
- Removes tree (or shrub) cohorts whose remaining density (resp. cover) is lower than a specified threshold.
- Merges surviving cohorts with recruitment in the
forestobject and prepares the input of functiongrowth()for the next annual time step. - Store current status of the
forestobject and update output tables/summaries.
18.4 Inputs and outputs
An important difference between fordyn() and the previous simulation functions is that it does not require a specific input object, as in spwb()or growth() functions. In other words, soil, vegetation, meteorology and control inputs are directly introduced as parameters to the function call to fordyn().
18.4.1 Soil, vegetation and meteorology
Soil
Soil input requirements are the same as for the former models and were fully described in section 2.3.
Vegetation
Unlike the former models, vegetation input for fordyn are objects of the class forest, which were described in section 2.4.4.
Metereological input
The minimum weather variables required to run the model are min/max temperatures (\(T_{min}\) and \(T_{max}\)), min/max relative humidity (\(RH_{min}\) and \(RH_{max}\)), precipitation (\(P\)) and solar radiation (\(Rad\)). Wind speed (\(u\)) is also needed, but the user may use missing values if not available (a default value will be used in this case). Wind speed is assumed to have been measured at a specific height above the canopy (by default at 2 m). Atmospheric \(CO_2\) concentration (\(C_{atm}\)) may also be specified, but if missing a default constant value is assumed, which is taken from the control parameters. Definitions and units of these variables were given in section 2.5.
18.4.2 Vegetation functional parameters
The forest dynamics model requires many functional parameters to be specified for plant cohorts. Some of them depend on whether the basic or advanced water balance is adopted, whereas others are inherited from the growth model. Here we report functional parameters needed in addition to those necessary for the growth model (see 15.5.2).
All of them concern the simulation of recruitment and are specified in the species parameter table (i.e. SpParams).
| Symbol | Units | R | Description |
|---|---|---|---|
| \(H_{seed}\) | \(cm\) | SeedProductionHeight |
Minimum height for seed production |
| \(SM\) | \(mg\) | SeedMass |
Seed dry mass |
| \(SL\) | \(yr\) | SeedLongevity |
Seedbank average longevity |
| \(Disp_{dist}\) | \(m\) | DispersalDistance |
Distance parameter for dispersal kernel |
| \(Disp_{shape}\) | DispersalShape |
Shape parameter for dispersal kernel | |
| \(P_{recr}\) | [0-1] | ProbRecr |
Probability of recruitment from seeds within the bioclimatic limits imposed by temperature, moisture and light thresholds |
| \(TCM_{recr}\) | \(^{\circ} \mathrm{C}\) | MinTempRecr |
Minimum average temperature (Celsius) of the coldest month for successful recruitment from seeds |
| \(MI_{recr}\) | MinMoistureRecr |
Minimum value of the moisture index (annual precipitation over annual PET) for successful recruitment from seeds | |
| \(FPAR_{recr}\) | % | MinFPARRecr |
Minimum percentage of PAR at the ground level for successful recruitment from seeds |
| \(N_{tree, recr}\) | \(ind \cdot ha^{-1}\) | RecrTreeDensity |
Density of tree recruits from seeds. |
| \(N_{tree, ingrowth}\) | \(ind \cdot ha^{-1}\) | IngrowthTreeDensity |
Density of trees reaching ingrowth DBH. |
| \(DBH_{tree,recr}\) | \(cm\) | RecrTreeDBH |
DBH for tree recruits from seeds or resprouting (e.g. 1 cm). |
| \(DBH_{tree,ingrowth}\) | \(cm\) | IngrowthTreeDBH |
Ingrowth DBH for trees (e.g. 7.5 cm). |
| \(H_{tree, recr}\) | \(cm\) | RecrTreeHeight |
Height for tree recruits from seeds or resprouting |
| \(Cover_{shrub, recr}\) | % | RecrShrubCover |
Recruitment cover for shrubs |
| \(H_{shrub, recr}\) | \(cm\) | RecrShrubHeight |
Recruitment height for shrubs |
| \(Z50_{recr}\) | mm | RecrZ50 |
Soil depth corresponding to 50% of fine roots for recruitment |
| \(Z95_{recr}\) | mm | RecrZ95 |
Soil depth corresponding to 95% of fine roots for recruitment |
| \(Resp_{fire}\) | RespFire |
Number of resprouts per stem after fire disturbance | |
| \(Resp_{dist}\) | RespDist |
Number of resprouts per stem after undefined disturbance (typically desiccation) | |
| \(Resp_{clip}\) | RespClip |
Number of resprouts per stem after clipping |
18.4.3 Control parameters
Control parameters modulate the overall behavior of fordyn simulations, which extend the parameters used for growth simulations (see section 15.5.3). First, there are parameters that regulate the application of seed production, seed bank dynamics, recruitment from seeds, resprouting and the removal of cohorts with few individuals:
applySeedBankDynamics [= TRUE]: Boolean flag to indicate that seed bank dynamics (seed production, seed bank mortality and seed rain) need to be simulated.applyRecruitment [= TRUE]: Boolean flag to indicate that recruitment from seeds is allowed.applyResprouting [= TRUE]: Boolean flag to indicate that resprouting is allowed.recruitmentMode [= "stochastic"]: String describing how recruitment from seeds is applied. Current accepted values are “deterministic” or “stochastic”.removeEmptyCohorts [= FALSE]: Boolean flag to indicate the removal of cohorts whose density is too low.minimumTreeCohortDensity [= 1]: Threshold of tree density resulting in cohort removal.minimumShrubCohortCover [= 0.01]: Threshold of shrub cover resulting in cohort removal.dynamicallyMergeCohorts [= TRUE]: Boolean flag to indicate that cohorts should be merged when possible. This option speeds up calculations but results in a loss of cohort identity and reinitialization of many state variables.
Next, a few parameters control the production/arrival of seeds:
seedRain [= NULL]: Vector of species names whose seed rain is to be added to seed bank, regardless of local seed production.seedProductionTreeHeight [= 300]: Default minimum tree height for producing seeds (when species parameterSeedProductionHeightis missing).seedProductionShrubHeight [= 30]: Default minimum shrub height for producing seeds (when species parameterSeedProductionHeightis missing).
Then we have default parameters determining whether recruitment occurs:
probRecr [= 0.05]: Default annual probability of recruitment (when species parameterProbRecris missing).minTempRecr [= 0]: Default threshold of minimum average temperature of the coldest month necessary for recruiting (when species parameterMinTempRecris missing).minMoistureRecr [= 0.3]: Default threshold of minimum moisture index (annual precipitation over annual ETP) necessary for recruiting (when species parameterMinMoistureRecris missing).minFPARRecr [= 10]: Default threshold of minimum fraction of PAR (in %) reaching the ground necessary for recruiting (when species parameterMinFPARRecris missing).
Finally, there are a set of parameters specifying default values for recruited cohort attributes:
recrTreeDBH [= 1]: Default DBH (cm) for recruited trees (when species parameterRecrTreeDBHis missing).recrTreeDensity [= 100]: Default density (ind·ha-1) for recruited trees (when species parameterRecrTreeDensityis missing).recrTreeHeight [= 100]: Default height (cm) for recruited trees (when species parameterRecrTreeHeightis missing).recrShrubCover [= 1]: Default cover (%) for recruited shrubs (when species parameterRecrShrubCoveris missing).recrShrubHeight [= 100]: Default height (cm) for recruited shrubs (when species parameterRecrShrubHeightis missing).recrTreeZ50 [= 100]: Default value for Z50 (mm) in recruited trees (when species parameterRecrZ50is missing).recrShrubZ50 [= 50]: Default value for Z50 (mm) in recruited shrubs (when species parameterRecrZ50is missing).recrTreeZ95 [= 1000]: Default value for Z95 (mm) in recruited trees (when species parameterRecrZ50is missing).recrShrubZ50 [= 500]: Default value for Z95 (mm) in recruited shrubs (when species parameterRecrZ50is missing).
18.4.4 Model output
| Element | Description |
|---|---|
StandSummary |
A data frame with stand-level summaries (leaf area index, tree basal area, tree density, shrub cover, etc.) at the beginning of the simulation and after each simulated year. |
SpeciesSummary |
A data frame with species-level summaries (leaf area index, tree basal area, tree density, shrub cover, etc.) at the beginning of the simulation and after each simulated year. |
CohortSummary |
A data frame with cohort-level summaries (leaf area index, tree basal area, tree density, shrub cover, etc.) at the beginning of the simulation and after each simulated year. |
TreeTable |
A data frame with tree-cohort data (species, density, diameter, height, etc.) at the beginning of the simulation (if any) and after each simulated year. |
DeadTreeTable |
A data frame with dead tree-cohort data (species, density, diameter, height, etc.) at the beginning of the simulation and after each simulated year. |
CutTreeTable |
A data frame with cut tree-cohort data (species, density, diameter, height, etc.) per each simulated year. |
ShrubTable |
A data frame with shrub-cohort data (species, density, cover, height, etc.) at the beginning of the simulation and after each simulated year. |
DeadShrubTable |
A data frame with dead shrub-cohort data (species, density, cover, height, etc.) at the beginning of the simulation (if any) and after each simulated year. |
CutShrubTable |
A data frame with cut shrub-cohort data (species, density, cover, height, etc.) per each simulated year. |
ForestStructures |
A list with the forest object of the stand at the beginning of the simulation and after each simulated year. |
GrowthResults |
A list with the results of calling function growth (i.e., see 15.5.4) for each simulated year. |
ManagementArgs |
If management is considered, a list of management arguments to be used in another call to fordyn(). |
NextInputObject |
An object of class growthInput to be used in a subsequent simulation. |
NextForestObject |
An object of class forest to be used in a subsequent simulation. |
18.5 Process details
18.5.1 Seed production and seed bank dynamics
The model considers mortality of seeds in the seed bank before adding new seeds. Annual seed bank mortality is simulated for each species using an exponential decay function, driven by the species-specific seed longevity (\(SL\)):
\[\begin{equation} Seeds_{t+1} = Seeds_{t} \cdot \exp(- 1/SL) \end{equation}\]
To determine seed production, the model determines the seed rain in the stand by determining which cohorts have total plant heights above maturity thresholds for trees and shrubs. This is done, for each cohort, by comparing its height with the species-specific parameter seed production height (\(H_{seed}\)). Missing values of \(H_{seed}\) are given values from control variables seedProductionTreeHeight and seedProductionShrubHeight for trees and shrubs, respectively. In addition to locally produced seeds, the user can use control parameter seedRain to specify a list of species names whose seeds arrive to the stand. Despite the origin, in non-spatial simulations the seed bank relative amount of species with seed rain is set to 100%, assuming that there is enough seed production to allow normal recruitment. However, when simulations are conducted in a spatial context, the species identity and relative amount of seeds is determined via dispersal sub-model (see chapter 20).
18.5.2 Recruitment from seeds
Actual sapling recruitment (i.e. recruitment of small trees, typically 1 cm diameter) depends on environmental conditions in the stand. Specifically, the model calculates, for the year that ended, the mean temperature of the coldest month (\(TCM\)), the moisture index (\(MI\)) and the fraction of photosynthetic active radiation reaching the ground, given the current structure (\(FPAR\)). These values are compared to species specific parameter thresholds \(TCM_{recr}\), \(MI_{recr}\) and \(FPAR_{recr}\), respectively. More specifically, a given species \(j\) can recruit only if \(TCM > TCM_{recr, j}\), \(MI > MI_{recr, j}\) and \(FPAR > FPAR_{recr, j}\). These filters do not ensure recruitment, as it is assumed that multiple processes can further determine the death of recruits. A probability of recruitment \(P_{recr}\) is used to represent these additional processes. When simulation of recruitment is stochastic, a given species will recruit if (in addition to the bioclimatic limits) a Bernouilli draw falls below \(P_{recr}\) and, if so, the initial density or cover will be fixed. If recruitment is deterministic, then the probability of recruitment is used to multiply the recruitment density or cover.
Tree recruitment density, diameter and height are determined by parameters \(N_{tree,recr}\), \(DBH_{tree, recr}\) and \(H_{tree, recr}\), respectively; whereas cover and height of shrub recruitment is determined by parameters \(Cover_{recr}\) and \(H_{shrub, recr}\), respectively. Density/cover of recruits is reduced depending on the relative amount of seeds in the seed bank (i.e. maximum density or cover will occur if seed bank levels are 100%). If stochastic simulation of recruitment is requested, then density or cover values are considered mean values of a Poisson distribution. Note that density of tree recruits will decrease during the years after recruitment due to a self-thinning process, until \(DBH\) reaches \(DBH_{tree,ingrowth}\) where \(N\) will be \(N_{tree,ingrowth}\), as explained in 17.4.1.
18.5.3 Resprouting
Currently, resprouting only occurs in the model if plant cohorts have been cut, burned or they died of desiccation (future model versions will include resprouting after fire impacts). In other words, resprouting does not occur for baseline mortality, self-thinning of recruits or due to starvation.
Resprouting survivorship
The model first determines resprouting survivorship, which depends on the disturbance type and, in the case of fire, on species identity. Unfortunately, we lack information on the mortality caused in resprouting by cutting or drought for the vast majority of species except for holm oak (Quercus ilex). Thus, we suggest applying tentatively to all species the mortality values obtained in different studies for this species: 2.5% after browsing (Espelta et al. 2006), 4% after cutting (Retana et al. 1992), 5% after drought (Lloret et al. 2004). Larger differences among species have been reported for fire (Espelta et al. 2012). This pattern may be linked to inter-specific differences in the size of the bud bank, the degree of bud protection or the amount of stored resources in belowground organs (burl, taproot). Accordingly, default values for survivorship are \(Resp_{clip} = 0.96\) (for cutting) and \(Resp_{dist} = 0.95\) (for dessication), whereas no default value is given for \(Resp_{fire}\). Naturally, the number of surviving trees (stumps) is given by the number of trees affected by the disturbance multiplied by the survivorship rate: \[\begin{equation} N_{surv} = N_{mort, dist} \cdot Resp_{dist} \end{equation}\] for a generic disturbance, but the same equation would apply for cutting, dessication or fire.
Resprouting vigor
For those individuals that survive the disturbance event and are able to resprout, the number of resprouts present in a given moment is a matter of the stump area and age. According to different literature sources, the number of resprouts initially produced is approximately 1.82 per 1 cm2 of stump area, but owing to the intra- and inter-individual competition (self-thinning), the number of resprouts will decrease with time following an exponential equation (Retana et al. 1999. Espelta et al. 1999, 2013):
\[\begin{equation} N_{resp} = N_{surv} \cdot \pi \cdot (DBH/2)^2 \cdot 1.82 \cdot 10^{-0.053\cdot x} \end{equation}\]
where \(DBH\) is diameter at breast height of the parent tree cohort, in cm, \(x = 5\) is time in years, \(N_{resp}\) is the number of resprouts, \(N_{surv}\) is the number of parent trees that survived.
Initial characteristics and self-thinning
Tree resprouts have all a starting diameter \(DBH_{tree, recr}\) and a height \(H_{tree, recr}\), but the root system is the same of the parent tree cohort. Analogously to recruitment from seeds, the initial density of tree resprout will decrease during the years after the initial resprouting due to a self-thinning process, until \(DBH\) reaches \(DBH_{tree,ingrowth}\) where \(N\) will be \(N_{tree,ingrowth}\), as explained in 17.4.1.
Resprouting in shrubs
In the case of shrub species, the same operations are done, but on cover values instead of density values. Shrub height of resprouts is assumed to be \(H_{shrub, recr}\).