The sitree package provides a framework to implement Single Tree forest growth models in a fast and memory efficient way. It keep tracks of all alive, dead, and removed trees in a robust and efficient way. SiTree is designed to run single tree simulations where trees can be defined by two time-dependent variables (such as diameter (or basal area), and height), and on time-independent variable, such as tree species. SiTree simulates birth, growth, and death of trees as well as management. Functions can also be defined that affect characteristics of the stand (external modifiers), such as climate change, or fertilization.
The easiest way to start with your own simulation is probably to modify the example functions provided (see the Test Equations vignette).
Two types of input are required by SiTree: tree
level and stand level. Tree level information is passed in
tree.df, while stand level information is passed in
stand.df.
tree.df should be a data frame with four columns named
plot.id, treeid, dbh,
height, and tree.sp, which correspond to a
stand/plot ID, a tree ID, diameter, height, and tree species.
Plot and stand data is passed in stand.df, which should
be a data frame or a list, with at least a column or element named
plot.id which should contain all the plot IDs present in
tree.df. Typical information provided in
stand.df are plot size, elevation, site index, plot
coordinates, distance to road, temperature or precipitation.
An example of tree data and stand data are provided.
## Loading required package: ggplot2## plot.id treeid dbh height tree.sp
## 1 91 108286 149 118 53
## 2 91 137120 133 107 53
## 3 91 108287 55 52 53
## 4 91 108268 235 151 49
## 5 91 137124 187 137 53
## 6 91 108271 162 124 49## $plot.id
## [1] 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
## [19] 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
## [37] 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54
## [55] 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
## [73] 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90
## [91] 91 92 93 94 95 96 97 98 99 100
##
## $SI.m
## [1] 11 11 8 11 8 8 6 14 11 11 8 23 6 6 11 11 17 14 8 8 8 6 17 8 6
## [26] 6 20 11 11 8 11 11 8 11 17 8 8 8 14 11 14 14 11 17 11 11 26 14 11 14
## [51] 11 11 6 14 17 8 8 11 8 17 14 14 14 14 17 8 8 8 14 11 11 20 17 17 11
## [76] 14 11 8 11 8 8 8 8 17 8 11 8 6 11 6 8 8 8 23 17 11 8 14 8 8
##
## $SI.spp
## [1] 2 2 3 3 3 3 3 3 3 3 2 1 2 3 3 3 1 3 3 3 3 2 3 3 2 2 2 2 3 2 2 2 2 2 1 2 2
## [38] 3 3 3 2 3 3 3 3 3 1 3 2 2 3 3 3 3 1 3 2 2 3 3 3 3 2 3 1 3 3 3 2 2 2 1 2 1
## [75] 3 3 3 2 3 3 2 2 3 3 3 3 2 3 3 2 3 3 3 1 3 2 3 1 2 1
##
## $prop.plot
## [1] 5 5 0 0 0 0 5 0 5 0 0 0 0 7 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4 0 0 0 0 0 0 0
## [38] 0 0 0 0 0 0 0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 7 3 0 0 0 0 0 0 6
## [75] 4 5 0 0 0 0 0 0 0 0 0 0 0 7 0 0 0 0 0 5 7 0 0 0 0 0
##
## $ha2total
## [1] 450.6190 450.6190 901.2380 901.2380 901.2380 901.2380 450.6190
## [8] 901.2380 450.6190 901.2380 901.2380 901.2380 901.2380 630.8666
## [15] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [22] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [29] 901.2380 360.4952 901.2380 901.2380 901.2380 901.2380 901.2380
## [36] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [43] 901.2380 901.2380 901.2380 270.3714 901.2380 901.2380 901.2380
## [50] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [57] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [64] 901.2380 901.2380 630.8666 270.3714 901.2380 901.2380 901.2380
## [71] 901.2380 901.2380 901.2380 540.7428 360.4952 450.6190 901.2380
## [78] 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380 901.2380
## [85] 901.2380 901.2380 901.2380 630.8666 901.2380 901.2380 901.2380
## [92] 2693.6380 901.2380 450.6190 630.8666 901.2380 901.2380 901.2380
## [99] 901.2380 901.2380
##
## $tree2ha
## [1] 80.00000 80.00000 40.00000 40.00000 40.00000 40.00000 80.00000
## [8] 40.00000 80.00000 40.00000 40.00000 40.00000 40.00000 57.14286
## [15] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [22] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [29] 40.00000 100.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [36] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [43] 40.00000 40.00000 40.00000 133.33333 40.00000 40.00000 40.00000
## [50] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [57] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [64] 40.00000 40.00000 57.14286 133.33333 40.00000 40.00000 40.00000
## [71] 40.00000 40.00000 40.00000 66.66667 100.00000 80.00000 40.00000
## [78] 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000 40.00000
## [85] 40.00000 40.00000 40.00000 57.14286 40.00000 40.00000 40.00000
## [92] 40.00000 40.00000 80.00000 57.14286 40.00000 40.00000 40.00000
## [99] 40.00000 40.00000Many of the main functions in a simulation use plot-level variables,
like competition indices (e.g. plot-level basal area). In order to make
the code more transparent, compact, and robust all variables required in
the sub-models that can be estimated from tree, stand and plot variables
are calculated in one place. The fn.prep.common.vars
function.
For example, if diameter increment is estimated as a function of
initial diameter, stand basal area and number of trees per ha, both
stand basal area variable and number of trees per ha should be
calculated in fn.prep.common.vars. Other typical examples
of variables calculated in the fn.prep.common.vars are top
height (the mean height of the trees with the largest diameter in a
stand), basal area of larger trees, or tree volume.
The fn.prep.common.vars function should be provided by
the user, so it fits the particular needs of the growth model selected
and the data. Calculating the most common variables used in forestry
should be straight forward as they are already provided in either the
SiTree or the SiTreeE package.
An example of a fn.prep.common.vars is provided in
SiTree.
## function (tr, fl, i.period, this.period, common.vars, vars.required,
## period.length, n.periods, ...)
## {
## if (length(common.vars) > 1)
## res <- common.vars
## else res <- list()
## others <- list(...)
## all.plot.vars <- data.table(plot.id = fl$plot.id, SI.m = fl$SI.m,
## SI.spp = fl$SI.spp, tree2ha = fl$tree2ha, ha2total = fl$ha2total,
## kom = fl$kom)
## res$i.stand <- match(tr$data[["plot.id"]], fl[["plot.id"]])
## res$i.tree <- match(fl$plot.id, tr$data$plot.id)
## res$tree.BA.m2 <- pi * (tr$data[["dbh.mm"]][, this.period]/1000/2)^2
## res$SBA.m2.ha <- tapply(res$tree.BA.m2 * fl[["tree2ha"]][res$i.stand],
## list(plot.id = tr$data[["plot.id"]]), FUN = sum)
## i.match.tapply <- match(tr$data[["plot.id"]], names(res$SBA.m2.ha))
## res$SBA.m2.ha <- as.vector(res$SBA.m2.ha[i.match.tapply])
## res$spp <- sp.classification(tree.sp = tr$data[["tree.sp"]],
## species.spruce = c(1, 2, 3), species.pine = c(10, 11,
## 20, 21, 29), species.harw = c(30, 31))
## all.tree.vars <- data.table(treeid = tr$data$treeid, plot.id = tr$data$plot.id,
## dbh.mm = tr$data$dbh.mm[, this.period], height.dm = tr$data$height.dm[,
## this.period], tree.sp = tr$data$tree.sp)
## all.tree.vars[all.plot.vars, `:=`(tree2ha, i.tree2ha), on = list(plot.id)]
## all.tree.vars[, `:=`(BA.m2, pi * (dbh.mm/1000/2)^2)]
## QMD.cm <- tapply(tr$data[["dbh.mm"]][, this.period], list(tr$data[["plot.id"]]),
## function(x.mm) {
## x.mm <- x.mm[is.finite(x.mm)]
## (sqrt(sum((x.mm/10)^2)/length(x.mm)))
## })
## res$QMD.cm <- as.vector(QMD.cm[i.match.tapply])
## tph <- tapply(fl[["tree2ha"]][res$i.stand], tr$data[["plot.id"]],
## sum)
## res$tph <- as.vector(tph[i.match.tapply])
## res$SDI <- res$tph * (res$QMD.cm/(10 * 2.54))^1.605
## pr.spp.ba <- data.frame(spru = rep(0, length(res$i.stand)),
## pine = 0, harw = 0, birch = 0, other = 0)
## pr.spp.ba$spru[res$spp == "spruce"] <- 1
## pr.spp.ba$pine[res$spp == "pine"] <- 1
## pr.spp.ba$birch[res$spp %in% c("birch")] <- 1
## pr.spp.ba$other[res$spp %in% c("other")] <- 1
## pr.spp.ba$harw[res$spp %in% c("birch", "other")] <- 1
## pr.spp.ba <- pr.spp.ba * res$tree.BA.m2
## dum.s <- tapply(pr.spp.ba$spru, tr$data$plot.id, sum)
## dum.p <- tapply(pr.spp.ba$pine, tr$data$plot.id, sum)
## dum.h <- tapply(pr.spp.ba$harw, tr$data$plot.id, sum)
## dum.b <- tapply(pr.spp.ba$birch, tr$data$plot.id, sum)
## dum.o <- tapply(pr.spp.ba$other, tr$data$plot.id, sum)
## pr.spp.ba <- data.frame(spru = as.vector(dum.s), pine = as.vector(dum.p),
## harw = as.vector(dum.h), birch = as.vector(dum.b), other = as.vector(dum.o))
## pr.spp.ba <- pr.spp.ba/with(pr.spp.ba, spru + pine + harw)
## pr.spp.ba <- pr.spp.ba[match(tr$data$plot.id, names(dum.s)),
## ] * 100
## res$pr.spp.ba <- pr.spp.ba
## rm(pr.spp.ba)
## res$PBAL.m2.ha <- ave(res$tree.BA.m2 * 10000/fl[["plot.size.m2"]][res$i.stand],
## tr$data$plot.id, FUN = function(X) {
## ord.x <- order(X)
## X <- sum(X[ord.x]) - cumsum(X[ord.x])
## X <- X[match(1:length(X), ord.x)]
## return(X)
## })
## previous.period <- paste0("t", i.period - 1)
## if ("stand.age.years" %in% names(fl)) {
## if (i.period == 0 & !is.data.frame(fl$stand.age.years)) {
## my.age <- fl$stand.age.years
## fl$stand.age.years <- data.frame(matrix(NA, ncol = n.periods,
## nrow = length(fl$plot.id)))
## names(fl$management) <- paste0("t", 1:n.periods)
## fl$stand.age.years[, this.period] <- my.age
## }
## if (i.period > 0 & "stand.age.years" %in% names(fl)) {
## fl$stand.age.years[, this.period] <- fl$stand.age.years[,
## previous.period] + 5
## stand.age.dt <- data.table(plot.id = fl$plot.id,
## SI.spp = fl$SI.spp, SI.m = fl$SI.m, stand.age.years = fl$stand.age.years[,
## this.period], waiting.time = 15)
## if (any(!is.na(fl$management[, this.period]))) {
## stand.age.dt[, `:=`(stands.ff, !substr(fl$management[,
## this.period], 1, 1) %in% c("0", "3"))]
## stand.age.dt[stands.ff == TRUE, `:=`(stand.age.years,
## -waiting.time + period.length/2)]
## stand.age.dt[, `:=`(dev.class, calculate.development.class(SI.spp = SI.spp,
## SI.m = SI.m, stand.age.years = stand.age.years))]
## }
## fl$stand.age.years[, this.period] <- stand.age.dt[,
## stand.age.years]
## if (any(is.na(stand.age.dt$stand.age.years)))
## browser()
## }
## res$dev.class <- calculate.development.class(SI.spp = fl$SI.spp,
## SI.m = fl$SI.m, stand.age.years = fl$stand.age.years[,
## this.period])
## res$tree.age <- data.table(age.years = fl$stand.age.years[,
## this.period][match(tr$data$plot.id, fl$plot.id)],
## treeid = tr$data$treeid)
## table(fl$stand.age.years[, this.period], useNA = "always")
## }
## vuprha.m3.ha <- NULL
## all.tree.vars[all.plot.vars, `:=`(kom, kom), on = "plot.id"]
## all.tree.vars[tree.sp == 12, `:=`(tree.sp, 10)]
## all.tree.vars[dbh.mm > 0, `:=`(c("vol.w.tr.m3", "vol.wo.tr.m3"),
## volume.norway(dbh.mm, height.dm, as.numeric(levels(tree.sp))[tree.sp],
## kom))]
## vuprha.m3.ha <- NULL
## all.tree.vars[all.plot.vars, `:=`(kom, i.kom), on = "plot.id"]
## all.tree.vars[tree.sp == 12, `:=`(tree.sp, 10)]
## all.tree.vars[dbh.mm > 0, `:=`(c("vol.w.tr.m3", "vol.wo.tr.m3"),
## volume.norway(dbh.mm, height.dm, as.numeric(levels(tree.sp))[tree.sp],
## kom))]
## all.tree.vars[, `:=`(vuprha.m3.ha, vol.wo.tr.m3 * tree2ha)]
## vuprha.m3.ha <- all.tree.vars[, sum(vuprha.m3.ha, na.rm = TRUE),
## by = plot.id]
## all.plot.vars[vuprha.m3.ha, `:=`(vuprha.m3.ha, V1), on = "plot.id"]
## all.plot.vars[is.na(vuprha.m3.ha), `:=`(vuprha.m3.ha, 0)]
## res$vuprha.m3.ha <- all.plot.vars$vuprha.m3.ha
## res$vol.wo.tr.m3.ha <- all.tree.vars$vuprha.m3.ha
## if (i.period == 0) {
## time.intern <- rep(NA, length.out = length(fl$plot.id))
## }
## else {
## time.intern <- fl$time.since.final.felling
## harv.last <- substr(fl$management[, this.period], 1,
## 1) %in% c("1", "2")
## time.intern[harv.last] <- period.length/2
## time.intern[!harv.last & !is.na(time.intern)] <- time.intern[!harv.last &
## !is.na(time.intern)] + period.length
## }
## fl$time.since.final.felling <- time.intern
## invisible(list(res = res, fl = fl))
## }
## <environment: namespace:sitree>trList and trListDead
classesIn order to efficiently store the list of all individual live and dead (and removed) trees, two Reference classes (or refclasses) are defined in SiTree. Refclases is chosen instead of S3 or S4 classes because refclasses objects are mutable and the usual R copy on modify semantics do not apply. When simulating for long periods, or for large datasets (e.g. a whole national forest inventory) the risk of running out of memory is not negligible. Using refclasses aim at maintaining the memory needs to the minimum by using mutable objects for storing the larger objects such as represented by the tree lists.
There are two Reference Classes implemented in the sitree package, one for live trees (trList) and other for dead trees (trListDead).
trList This class has two fields, data and nperiods. Under data basic information for each tree is stored (a unique stand ID plot.id, a unique tree ID treeid, and dbh and height for each period, dbh.mm, and height.dm). The nperiods field is an integer that stores the number of periods to be simulated.
trListDead This class extends trList. DBH and heights measured while the tree was alive are stored under the data field. Also under this field information on how long the tree has been in the simulation can be found (yrs.sim). In this class the new field last.measurement stores the dimensions of the tree when it died or was removed. How these dimensions are calculated is defined on the dead.trees.growth function.
Reference Classes objects are mutable, they don’t use R’s usual copy-on-modify semantics, but are modified in place.
We have provided a function to convert the sitree()
output containing trList and *trListDead** objects to a data
frame, the sitree2dataframe function. The resulting data
frame follows the usual R copy on modify semantics, and most users might
be more comfortable with it.
result.sitree <- sitree (tree.df = stand.west.tr,
stand.df = stand.west.st,
functions = list(
fn.growth = 'grow.dbhinc.hgtinc',
fn.mort = 'mort.B2007',
fn.recr = 'recr.BBG2008',
fn.management = NULL,
fn.tree.removal = NULL,
fn.modif = NULL,
fn.prep.common.vars = 'prep.common.vars.fun'
),
n.periods = 12,
period.length = 5,
mng.options = NA,
print.comments = FALSE,
fn.dbh.inc = 'dbhi.BN2009',
fn.hgt.inc = 'height.korf'
)
str(result.sitree$live)## Reference class 'trList' [package "sitree"] with 2 fields
## $ data :List of 6
## ..$ plot.id : num [1:879] 2 2 2 2 2 2 2 2 2 1 ...
## ..$ treeid : int [1:879] 4 5 6 7 8 9 10 11 13 16 ...
## ..$ dbh.mm :'data.frame': 879 obs. of 13 variables:
## .. ..$ t0 : num [1:879] 50 50 116 50 52 84 68 90 115 116 ...
## .. ..$ t1 : num [1:879] 61 61 138 61 64 101 83 108 137 137 ...
## .. ..$ t2 : num [1:879] 71 71 157 71 75 116 96 124 156 155 ...
## .. ..$ t3 : num [1:879] 80 80 174 80 85 130 108 139 173 171 ...
## .. ..$ t4 : num [1:879] 89 89 189 89 94 143 119 152 188 186 ...
## .. ..$ t5 : num [1:879] 97 97 203 97 103 155 129 164 202 200 ...
## .. ..$ t6 : num [1:879] 105 105 216 105 111 166 139 176 215 213 ...
## .. ..$ t7 : num [1:879] 113 113 229 113 119 177 149 187 228 225 ...
## .. ..$ t8 : num [1:879] 121 121 241 121 127 187 158 198 240 236 ...
## .. ..$ t9 : num [1:879] 129 129 252 129 135 197 167 208 251 247 ...
## .. ..$ t10: num [1:879] 136 136 263 136 143 207 176 218 262 257 ...
## .. ..$ t11: num [1:879] 143 143 273 143 150 216 184 228 272 267 ...
## .. ..$ t12: num [1:879] 150 150 283 150 157 225 192 237 282 277 ...
## ..$ height.dm:'data.frame': 879 obs. of 13 variables:
## .. ..$ t0 : num [1:879] 8 42 113 53 53 81 81 104 113 113 ...
## .. ..$ t1 : num [1:879] 18 52 128 63 64 94 94 118 128 127 ...
## .. ..$ t2 : num [1:879] 27 61 140 72 74 105 104 129 140 138 ...
## .. ..$ t3 : num [1:879] 35 69 150 80 82 115 113 139 150 147 ...
## .. ..$ t4 : num [1:879] 42 76 158 87 89 123 121 147 158 155 ...
## .. ..$ t5 : num [1:879] 48 82 165 93 96 130 128 154 165 162 ...
## .. ..$ t6 : num [1:879] 54 88 171 99 102 136 135 161 171 168 ...
## .. ..$ t7 : num [1:879] 60 94 177 105 108 142 141 167 177 174 ...
## .. ..$ t8 : num [1:879] 66 100 182 111 113 147 146 173 182 179 ...
## .. ..$ t9 : num [1:879] 71 105 187 116 118 152 151 178 187 184 ...
## .. ..$ t10: num [1:879] 76 110 192 121 123 157 156 183 192 188 ...
## .. ..$ t11: num [1:879] 80 114 196 125 127 161 160 188 196 192 ...
## .. ..$ t12: num [1:879] 84 118 200 129 131 165 164 192 200 196 ...
## ..$ yrs.sim : num [1:879] 60 60 60 60 60 60 60 60 60 60 ...
## ..$ tree.sp : Factor w/ 29 levels "1","2","3","10",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ nperiods: int 12
## and 19 methods, of which 5 are possibly relevant:
## addTrees, as.list, extractTrees, getTrees, show#envRefClass## treeid plot.id tree.sp period dbh.mm height.dm
## 1 4 2 1 t0 50 8
## 2 5 2 1 t0 50 42
## 3 6 2 1 t0 116 113
## 4 7 2 1 t0 50 53
## 5 8 2 1 t0 52 53
## 6 9 2 1 t0 84 81The sitree() function is the core function of the
SiTree package. It is the function that runs the
simulations. It requires tree data (tree.df), stand/plot
data (stand.df), and a list of functions to be used in the
simulation (functions), the number of periods for which to
run the simulation (n.periods), and the period length
(period.length). Management options can be passed through
the mng.options argument, and it is also possible to print
comments about the progress of the simulation selecting
print.comments = TRUE. Additional arguments needed by the
selected functions go into the ellipsis (‘…’) and can be retrieved by
simply converting it to a list,
e.g. arguments <- list(...).
The functions argument must be a list containing at
least 7 elements:
fn.growth the name of the function that implements the
growth sub-modelfn.mort the name of the function that implements the
mortality sub-modelfn.recr the name of the function that implements the
recruitment sub-model: ingrowth, and natural and artificial
regenerationfn.management the name of the function that implements
the management (e.g. harvest), can be NULLfn.tree.removal the name of the function that
implements the management at the tree level, that is the selection of
trees to be fell, can be NULL.fn.modif the name of the function to calculate the
effect of external modifiers, can be NULLfn.prep.common.vars the name of the function to
calculate auxiliary variables, such as basal area of the stand (see ‘The
fn.prep.common.vars’ subsection above) .Further details on the requirements of the functions listed above can be found under the section “The user-defined functions”.
The sitree() function is a flexible framework for forest
growth simulations. Any growth sub-model, mortality sub-model,
management, etc. can be used. Some examples are provided in
SiTree and in SiTreeE, but generally,
the submodels functions need to be provided by the user. The examples
provided in SiTree and in SiTreeE can
be used as a template. To debug the user-defined functions we suggest to
use the provided example as a starting point, set
print.comments = TRUE and switch the submodels functions
one by one to test them.
An example of how the list provided in the functions
argument of sitree should look like is given below, and
further details on each of the functions are provided next.
The fn.growth function should return a data frame with
two columns giving diameter increment (dbh.inc.mm) and
height increment (hgt.inc.dm) of all live trees. This data
frame should only contain numerical data (no missing data allowed). Care
must be taken to ensure that the order matches that of the tree list.
Examples of the growth functions are provided as
grow.dbhinc.hgtinc, dbhi.BN2009, and
height.korf.
The fn.mort function should return a
TRUE/FALSE vector of same length as the number
of trees in the tree list. TRUE indicates a tree that will
die before the next period, and FALSE indicates a tree that
will stay alive. An example of a fn.mort function is
provided in mort.B2007.
The fn.recr function is the function that estimates
recruitment, the new trees for the next period. This function should
return a list of new trees (or an empty list if there are no new trees)
with elements plot.id, treeid,
dbh.mm, height.dm, yrs.sim
(indicates when are the trees incorporated to the plot, for example, in
the middle of the period), and tree.sp.An example of a
fn.recr function is provided in
recr.BBG2008.
fn.management is optional. It should return a list, with
at least one element called management which should be a
vector with length equal to the number of plots in
stand.df. The example we provide uses a simple code to
define management (a five characters string indicating with a binary
code (1 = present, 0 = absent) the treatments to be executed:
harvest-thinning-fertilization-pruning-other), but any other way to code
management can be used, as far as fn.management returns a
vector. When no management will take place during the simulation
fn.management can be set to NULL. An example
of a fn.management function is provided in
management.prob.
fn.tree.removal is optional. It should return a
TRUE/FALSE vector indicating which trees are
to be removed. The vector should have the same length as the number of
trees alive at the current period. When no tree removal will take place
during the simulation (no harvest is allowed)
fn.tree.removal can be set to NULL. An example of a
fn.tree.removal function is provided in
mng.tree.removal.
fn.modif is a function that can be used to modify
characteristics of the plot or stand, such as site index. This function
is optional, and no example is provided in the current version of the
package. It should return a list with names matching some of those in
the stand.df data frame. After the external modifiers are
calculated with the function defined as fn.modif, the
elements in the plot data that matches those of the results of
fn.modif are replaced before the rest of the simulation
continues. For example, if the plot has been fertilized and we can
assume that SI has increased by 2 meters, the fn.modif
function needs to return a list with a SI element with all SI as in the
plot data frame except for those that have changed.
fn.prep.common.vars is the function used to calculate
everything needed for the fn.growth, fn.mort,
etc to be calculated. For example, the fn.prep.common.vars
function is the place to calculate stand competition indices, volume,
stand age, etc. An example is given in the function
prep.common.vars.fn.