sensRange {FME}  R Documentation 
Sensitivity Ranges of a Timeseries or 1D Variables
Description
Given a model consisting of differential equations, estimates the global effect of certain (sensitivity) parameters on a time series or on 1D spatial series of selected sensitivity variables.
This is done by drawing parameter values according to some predefined distribution, running the model with each of these parameter combinations, and calculating the values of the selected output variables at each output interval.
This function thus produces 'envelopes' around the sensitivity variables.
Usage
sensRange(func, parms = NULL, sensvar = NULL, dist = "unif",
parInput = NULL, parRange = NULL, parMean = NULL,
parCovar = NULL, map = 1, num = 100, ...)
## S3 method for class 'sensRange'
summary(object, ...)
## S3 method for class 'summary.sensRange'
plot(x, xyswap = FALSE,
which = NULL, legpos = "topleft",
col = c(grey(0.8), grey(0.7)),
quant = FALSE, ask = NULL, obs = NULL,
obspar = list(), ...)
## S3 method for class 'sensRange'
plot(x, xyswap = FALSE,
which = NULL, ask = NULL, ...)
Arguments
func 
an Rfunction that has as first argument 
parms 
parameters passed to 
sensvar 
the output variables for which the sensitivity needs
to be estimated. Either 
dist 
the distribution according to which the parameters should be
generated, one of 
parRange 
the range (min, max) of the sensitivity parameters, a
matrix or (preferred) a data.frame with one row for each parameter,
and two columns with the minimum (1st) and maximum (2nd) value. The
rownames of 
parInput 
a matrix with dimension (*, npar) with the values of the sensitivity parameters. 
parMean 
only when 
parCovar 
only when 
num 
the number of times the model has to be run. Set large
enough. If 
map 
the column number with the (independent) mapping variable
in the output matrix returned by 
object 
an object of class 
x 
an object of class 
legpos 
position of the legend; set to 
xyswap 
if 
which 
the name or the index to the variables that should be plotted. Default = all variables. 
col 
the two colors of the polygons that should be plotted. 
quant 
if 
ask 
logical; if 
obs 
a 
obspar 
additional graphics arguments passed to 
... 
additional arguments passed to 
Details
Models solved by integration (i.e. by using one of 'ode', 'ode.1D',
'ode.band', 'ode.2D'
), have the output already in a form usable by
sensRange
.
Value
a data.frame
of type sensRange
containing the parameter
set and the corresponding values of the sensitivity output variables.
The list returned by sensRange
has a method for the generic
functions summary
,plot
and
plot.summary
– see note.
Note
The following methods are included:

summary, estimates summary statistics for the sensitivity variables, a data.frame with as many rows as there are mapping variables (or rows in the matrix returned by
func
) and the following columns:x
, the mapping value,Mean
, the mean,sd
, the standard deviation,Min
, the minimal value,Max
, the maximal value,q25
,q50
,q75
, the 25th, 50 and 75% quantile 
plot, produces a "matplot" of the
sensRange
output, one plot for each sensitivity variable and with the mapping variable on the xaxis.Each variable will be plotted in a separate figure, and the figures aligned in a rectangular grid, unless par
mfrow
is passed as an argument. summary.plot, produces a plot of the summary of the
sensRange
output, one plot for each sensitivity variable and with the ranges and mean + standard deviation or the quantiles as coloured polygons.Each variable will be plotted in a separate figure, and the figures aligned in a rectangular grid, unless par
mfrow
is passed as an argument.
The output for models solved by a steadystate solver (i.e. one of
'steady', 'steady.1D', 'steady.band', 'steady.2D'
, needs to be
rearranged – see examples.
For plot.summary.sensRange
and plot.sensRange
,
the number of panels per page is automatically determined up to 3 x 3
(par(mfrow = c(3, 3))
). This default can be overwritten by
specifying userdefined settings for mfrow
or mfcol
.
Set mfrow
equal to NULL
to avoid the plotting function to
change userdefined mfrow
or mfcol
settings.
Other graphical parameters can be passed as well. Parameters
are vectorized, either according to the number of plots
(xlab, ylab
, main, sub
, xlim, ylim
, log
,
asp, ann, axes, frame.plot
,panel.first,panel.last
,
cex.lab,cex.axis,cex.main
) or
according to the number of lines within one plot (other parameters
e.g. col
, lty
, lwd
etc.) so it is possible to
assign specific axis labels to individual plots, resp. different plotting
style. Plotting parameter ylim
, or xlim
can also be a list
to assign different axis limits to individual plots.
Similarly, the graphical parameters for observed data, as passed by
obspar
can be vectorized, according to the number of observed
data sets (when obs
is a list
).
The data.frame
of type sensRange
has several attributes,
which remain hidden, and which are generally not of practical use
(they are needed for the S3 methods).
There is one exception, i.e. if parameter values are imposed via
argument parInput
, and these parameters are generated by a
Markov chain (modMCMC
). If the number of draws,
num
, is less than the number of rows in parInput
, then
num
random draws will be taken. Attribute, "pset" then contains
the index to the parameters that have been selected.
The sensRange
method only represents the distribution of the
model response variables as a function of the parameter values. But an
additional source of noise is due to the model error, as
represented by the sampled values of sigma in the Markov chain. In
order to represent also this source of error, gaussian noise should be
added to each sensitivity output variables, with a standard deviation
that corresponds to the original parameter draw – see vignette
"FMEother".
Author(s)
Karline Soetaert <karline.soetaert@nioz.nl>
References
Soetaert, K. and Petzoldt, T. 2010. Inverse Modelling, Sensitivity and Monte Carlo Analysis in R Using Package FME. Journal of Statistical Software 33(3) 1–28. doi:10.18637/jss.v033.i03
Examples
## =======================================================================
## Bacterial growth model from Soetaert and Herman, 2009
## =======================================================================
pars < list(gmax = 0.5,eff = 0.5,
ks = 0.5, rB = 0.01, dB = 0.01)
solveBact < function(pars) {
derivs < function(t,state,pars) { # returns rate of change
with (as.list(c(state,pars)), {
dBact < gmax*eff * Sub/(Sub + ks)*Bact  dB*Bact  rB*Bact
dSub < gmax * Sub/(Sub + ks)*Bact + dB*Bact
return(list(c(dBact,dSub)))
})
}
state < c(Bact = 0.1,Sub = 100)
tout < seq(0, 50, by = 0.5)
## ode solves the model by integration ...
return(as.data.frame(ode(y = state, times = tout, func = derivs,
parms = pars)))
}
out < solveBact(pars)
mf <par(mfrow = c(2,2))
plot(out$time, out$Bact, main = "Bacteria",
xlab = "time, hour", ylab = "molC/m3", type = "l", lwd = 2)
## the sensitivity parameters
parRanges < data.frame(min = c(0.4, 0.4, 0.0), max = c(0.6, 0.6, 0.02))
rownames(parRanges)< c("gmax", "eff", "rB")
parRanges
tout < 0:50
## sensitivity to rB; equallyspaced parameters ("grid")
SensR < sensRange(func = solveBact, parms = pars, dist = "grid",
sensvar = "Bact", parRange = parRanges[3,], num = 50)
Sens <summary(SensR)
plot(Sens, legpos = "topleft", xlab = "time, hour", ylab = "molC/m3",
main = "Sensitivity to rB", mfrow = NULL)
## sensitivity to all; latin hypercube
Sens2 < summary(sensRange(func = solveBact, parms = pars, dist = "latin",
sensvar = c("Bact", "Sub"), parRange = parRanges, num = 50))
## Plot all variables; plot mean + sd, min max
plot(Sens2, xlab = "time, hour", ylab = "molC/m3",
main = "Sensitivity to gmax,eff,rB", mfrow = NULL)
par(mfrow = mf)
## Select one variable for plotting; plot the quantiles
plot(Sens2, xlab = "time, hour", ylab = "molC/m3", which = "Bact", quant = TRUE)
## Add data
data < cbind(time = c(0,10,20,30), Bact = c(0,1,10,45))
plot(Sens2, xlab = "time, hour", ylab = "molC/m3", quant = TRUE,
obs = data, obspar = list(col = "darkblue", pch = 16, cex = 2))