Forest growth

About this vignette

This document describes how to run the growth model of medfate, implemented in function growth(). All the details of the model design and formulation can be found at https://emf-creaf.github.io/medfatebook/index.html.

Preparing model inputs

Model inputs are explained in greater detail in vignette ‘Simulation inputs’. Here we only review the different steps required to run function growth().

Soil, vegetation, meteorology and species data

Soil information needs to be entered as a data frame with soil layers in rows and physical attributes in columns. Soil physical attributes can be initialized to default values, for a given number of layers, using function defaultSoilParams():

spar = defaultSoilParams(2)

The soil input for water balance simulation is actually a list of class soil that is created using a function with the same name:

examplesoil = soil(spar)

As explained in the package overview, models included in medfate were primarily designed to be ran on forest inventory plots. Here we use the example object provided with the package:

data(exampleforestMED)
exampleforestMED

## $ID
## [1] "1"
## 
## $patchsize
## [1] 10000
## 
## $treeData
##   Species   N   DBH Height Z50  Z95
## 1     148 168 37.55    800 750 3000
## 2     168 384 14.60    660 750 3000
## 
## $shrubData
##   Species Cover Height Z50  Z95
## 1     165  3.75     80 300 1500
## 
## $herbCover
## [1] 10
## 
## $herbHeight
## [1] 20
## 
## attr(,"class")
## [1] "forest" "list"

In the basic water balance, only mean temperature, precipitation and potential evapotranspiration is required, but radiation may also be necessary to simulate snow melt.

data(examplemeteo)
head(examplemeteo)

##            MeanTemperature MinTemperature MaxTemperature Precipitation
## 2001-01-01      3.57668969     -0.5934215       6.287950      4.869109
## 2001-01-02      1.83695972     -2.3662458       4.569737      2.498292
## 2001-01-03      0.09462563     -3.8541036       2.661951      0.000000
## 2001-01-04      1.13866156     -1.8744860       3.097705      5.796973
## 2001-01-05      4.70578690      0.3288287       7.551532      1.884401
## 2001-01-06      4.57036721      0.5461322       7.186784     13.359801
##            MeanRelativeHumidity MinRelativeHumidity MaxRelativeHumidity
## 2001-01-01             78.73709            65.15411           100.00000
## 2001-01-02             69.70800            57.43761            94.71780
## 2001-01-03             70.69610            58.77432            94.66823
## 2001-01-04             76.89156            66.84256            95.80950
## 2001-01-05             76.67424            62.97656           100.00000
## 2001-01-06             89.01940            74.25754           100.00000
##            Radiation WindSpeed WindDirection       PET
## 2001-01-01  12.89251  2.000000           172 1.3212770
## 2001-01-02  13.03079  7.662544           278 2.2185985
## 2001-01-03  16.90722  2.000000           141 1.8045176
## 2001-01-04  11.07275  2.000000           172 0.9200627
## 2001-01-05  13.45205  7.581347           321 2.2914449
## 2001-01-06  12.84841  6.570501           141 1.7255058

Finally, simulations in medfate require a data frame with species parameter values, which we load using defaults for Catalonia (NE Spain):

data("SpParamsMED")

Simulation control

Apart from data inputs, the behaviour of simulation models can be controlled using a set of global parameters. The default parameterization is obtained using function defaultControl():

control = defaultControl("Granier")

Here we will run growth simulations using the basic water balance model (i.e. transpirationMode = "Granier"). The complexity of the soil water balance calculations can be changed by using "Sperry" as input to defaultControl().

Growth input object

A last object is needed before calling simulation functions, called growthInput. This is analogous to spwbInput and consists in the compilation of aboveground and belowground parameters and the specification of additional parameter values for each plant cohort. The object can be obtained by calling function growthInput(), but if one has a forest object, it can be generated more directly using function forest2growthInput():

x = forest2growthInput(exampleforestMED, examplesoil, SpParamsMED, control)

All the input information for forest data and species parameter values can be inspected by printing different elements of the input object, which are:

names(x)

##  [1] "control"             "soil"                "canopy"             
##  [4] "cohorts"             "above"               "below"              
##  [7] "belowLayers"         "paramsPhenology"     "paramsAnatomy"      
## [10] "paramsInterception"  "paramsTranspiration" "paramsWaterStorage" 
## [13] "paramsGrowth"        "paramsAllometries"   "internalPhenology"  
## [16] "internalWater"       "internalCarbon"      "internalAllocation" 
## [19] "internalMortality"   "internalRings"

First, information about the cohort species is found in element cohorts (i.e. code, species and name):

x$cohorts

##         SP              Name
## T1_148 148  Pinus halepensis
## T2_168 168      Quercus ilex
## S1_165 165 Quercus coccifera

Element above contains the aboveground structure data that we already know, but with additional columns that describe the estimated initial amount of sapwood area and carbon reserves (if required):

x$above

##         SP        N   DBH Cover   H        CR          SA   LAI_live
## T1_148 148 168.0000 37.55    NA 800 0.6605196 349.0550668 1.00643723
## T2_168 168 384.0000 14.60    NA 660 0.6055642  97.2338913 0.92661573
## S1_165 165 749.4923    NA  3.75  80 0.8032817   0.8561788 0.03965932
##        LAI_expanded LAI_dead
## T1_148   1.00643723        0
## T2_168   0.92661573        0
## S1_165   0.03965932        0

As with the soil input object, the growthInput object will be modified during simulations. In the case of growth(), these modifications concern SA, fastCstorage, slowCstorage, LAI_expanded and LAI_dead in element above, as well as canopy variables.

Elements starting with params* contain cohort-specific model parameters. Aboveground parameters related to plant transpiration can be seen in paramsTransp:

x$paramsTransp

##        Tmax_LAI Tmax_LAIsq Psi_Extract Psi_Critic      WUE pRootDisc
## T1_148    0.134     -0.006   -2.600000  -3.816852 3.613308         0
## T2_168    0.134     -0.006   -3.001667  -3.001667 1.475715         0
## S1_165    0.134     -0.006   -3.455000  -6.950920 2.433101         0

Belowground parameters can be seen in below and include root distribution as well as maximum root and rhizosphere conductances by soil layer:

x$below

##        Z50  Z95 coarseRootSoilVolume
## T1_148 750 3000            0.6815696
## T2_168 750 3000            0.6815696
## S1_165 300 1500            0.4242041

Parameters related to plant growth can be seen in paramsGrowth and paramsAllometries:

x$paramsGrowth

##           RERleaf   RERsapwood RERfineroot RGRleafmax RGRsapwoodmax
## T1_148 0.00260274 6.849315e-05 0.002054795       0.01         0.002
## T2_168 0.00260274 6.849315e-05 0.002054795       0.01         0.002
## S1_165 0.00260274 6.849315e-05 0.002054795       0.01         0.002
##        RGRfinerootmax    SRsapwood  SRfineroot fHDmin fHDmax     WoodC
## T1_148            0.1 0.0001261398 0.001897231     26    100 0.4978954
## T2_168            0.1 0.0001261398 0.001897231     16    160 0.4740096
## S1_165            0.1 0.0001261398 0.001897231     NA     NA 0.4743383

x$paramsAllometries

##             Aash     Bash      Absh      Bbsh     Acr   B1cr     B2cr
## T1_148 5.8458000 1.494400 0.3596000 0.7138000 1.99500 -0.649 -0.02000
## T2_168 1.8574862 1.885548 0.5238830 0.7337293 1.98539 -0.552 -0.01386
## S1_165 0.1305509 2.408443 0.5147731 0.5311554 1.50600 -0.706 -0.07800
##                B3cr     C1cr     C2cr       Acw    Bcw
## T1_148 -0.000120000 -0.00400 -0.15900 0.6415296 0.7310
## T2_168 -0.000110736 -0.00685 -0.20101 0.5681897 0.7974
## S1_165  0.000180000 -0.00700  0.00000 0.8390000 0.7350

Finally, elements starting with internal* contain state variables required to keep track of plant status. For example, the carbon storage levels can be seen in internalCarbon:

x$internalCarbon

##        sugarLeaf  starchLeaf sugarSapwood starchSapwood
## T1_148 0.4029239 0.002567561    0.5717221    0.01211977
## T2_168 0.3585751 0.002567561    1.0288757    0.01198800
## S1_165 0.7223526 0.002567561   -0.1709299    0.01075605

Executing the growth model

In this vignette we will use one year of example meteorological data with the following precipitation and temperature seasonal trends:

G1<-growth(x, examplemeteo, latitude = 41.82592, elevation = 100)

## Initial soil water content (mm): 185.069
## Initial snowpack content (mm): 0
## Performing daily simulations 
##  Year 2001:.................................... [update structural variables] done.
## Final soil water content (mm): 166.178
## Final snowpack content (mm): 0
## Change in soil water content (mm): -18.891
## Soil water balance result (mm): -18.891
## Change in snowpack water content (mm): 0
## Snowpack water balance result (mm): 7.10543e-15
## Water balance components:
##   Precipitation (mm) 513
##   Rain (mm) 462 Snow (mm) 51
##   Interception (mm) 92 Net rainfall (mm) 371
##   Infiltration (mm) 412 Runoff (mm) 10 Deep drainage (mm) 101
##   Soil evaporation (mm) 41 Transpiration (mm) 289

While water/carbon balances and growth of leaf/sapwood/root tissues is performed on a daily basis, structural variables are only updated once a year, here coinciding with the end of the simulation.

Plotting results

In addition to all the plots that were available for result of simulations with function spwb(), one can display the following information out of the output of function growth():

Carbon balance components

To inspect components of the plant carbon balance we can first display gross photosynthesis expressed as assimilation relative to biomass:

plot(G1, "GrossPhotosynthesis", bySpecies = T)

Then we can draw the respiration costs (maintenance and growth) in the same units:

plot(G1, "MaintenanceRespiration", bySpecies = T)

plot(G1, "GrowthCosts", bySpecies = T)

Finally we can display the daily negative or positive balance of the plant, which determines changes in plant carbon pools:

plot(G1, "CarbonBalance", bySpecies = T)

Carbon storage

Carbon assimilation and respiration rates define the dynamics of stored carbon, which in medfate is divided into two organs (leaves and sapwood) and two chemical compounds (soluble sugars and starch):

plot(G1, "SugarLeaf", bySpecies = T)

plot(G1, "StarchLeaf", bySpecies = T)

plot(G1, "SugarSapwood", bySpecies = T)

plot(G1, "StarchSapwood", bySpecies = T)

Leaf and sapwood area changes

Leaf and sapwood area dynamics arising from the interplay between growth and senescence of tissues can be inspected using:

plot(G1, "LeafArea", bySpecies = T)

plot(G1, "SapwoodArea", bySpecies = T)

The corresponding daily growth rates can also be displayed:

plot(G1, "LAgrowth", bySpecies = T)

plot(G1, "SAgrowth", bySpecies = T)

Finally, recall that one can interactively create plots using function shinyplot, e.g.:

shinyplot(G1)