plot.ODEmorris {ODEsensitivity}

Plot of the Results of Morris Screening for Objects of Class ODEmorris

Description

plot.ODEmorris plots the results of Morris screening for objects of class ODEmorris.

Usage

## S3 method for class 'ODEmorris'
plot(x, pars_plot = NULL, state_plot = names(x)[1],
  kind = "sep", colors_pars = NULL, main_title = NULL,
  legendPos = "outside", type = "l", ...)

Arguments

x

[ODEmorris]
output of ODEmorris (of class ODEmorris).

pars_plot

[character(k)]
names of the k parameters to be plotted. If NULL (the default), all parameters are plotted.

state_plot

[character(1)]
name of the state variable to be plotted. Defaults to the name of the first state variable.

kind

[character(1)]
kind of the plot, choose between "sep" and "trajec" (see details).

colors_pars

[character(>= k)]
vector of the colors to be used for the k different parameters. Must be at least of length k (only the first k elements will be used, though). If NULL (the default), rainbow(k) is used.

main_title

[character(1)]
title for the plot. If kind = "sep", this is the overall title for the two separate plots. If NULL (the default), a standard title is generated.

legendPos

[character(1)]
keyword for the legend position, either one of those specified in legend or "outside" (the default), which means the legend is placed under the plot (useful, if there are many parameters in the model).

type

[character(1)]
plot type, i.e. "p", "l", "b", "c", "o", "s", "h" or "n". Defaults to "l".

...

additional arguments passed to plot.default.

Details

Morris sensitivity indices are plotted for one state variable (chosen by argument state_plot) and the parameters named in pars_plot. If no parameters are named in pars_plot, the sensitivity indices for all parameters are plotted. There are two kinds of plots:

• kind = "sep": separate plots of the Morris sensitivity indices \mu^* and \sigma against time

• kind = "trajec": plot of \mu^* against \sigma

Value

TRUE (invisible; for testing purposes).

Note

Not all plotting arguments can be passed by ..., for example xlab and ylab are fixed.

Author(s)

Stefan Theers, Frank Weber

See Also

ODEmorris, morris

Examples

##### Lotka-Volterra equations #####
LVmod <- function(Time, State, Pars) {
  with(as.list(c(State, Pars)), {
    Ingestion    <- rIng  * Prey * Predator
    GrowthPrey   <- rGrow * Prey * (1 - Prey/K)
    MortPredator <- rMort * Predator
    
    dPrey        <- GrowthPrey - Ingestion
    dPredator    <- Ingestion * assEff - MortPredator
    
    return(list(c(dPrey, dPredator)))
  })
}
LVpars  <- c("rIng", "rGrow", "rMort", "assEff", "K")
LVbinf <- c(0.05, 0.05, 0.05, 0.05, 1)
LVbsup <- c(1.00, 3.00, 0.95, 0.95, 20)
LVinit  <- c(Prey = 1, Predator = 2)
LVtimes <- c(0.01, seq(1, 50, by = 1))
set.seed(7292)
# Warning: The following code might take very long!

LVres_morris <- ODEmorris(mod = LVmod,
                          pars = LVpars,
                          state_init = LVinit,
                          times = LVtimes,
                          binf = LVbinf,
                          bsup = LVbsup,
                          r = 500,
                          design = list(type = "oat", 
                                        levels = 10, grid.jump = 1),
                          scale = TRUE,
                          ode_method = "lsoda",
                          parallel_eval = TRUE,
                          parallel_eval_ncores = 2)
my_cols <- c("firebrick", "orange2", "dodgerblue", 
             "forestgreen", "black")
plot(LVres_morris, kind = "sep", colors_pars = my_cols)
plot(LVres_morris, pars_plot = c("rGrow", "rMort"), state_plot = "Predator", 
     kind = "trajec", colors_pars = my_cols[2:3])


##### A network of 4 mechanical oscillators connected in a circle #####
M_mat <- rep(2, 4)
K_mat <- diag(rep(2 * (2*pi*0.17)^2, 4))
K_mat[1, 2] <- K_mat[2, 3] <- 
  K_mat[3, 4] <- K_mat[1, 4] <- 2 * (2*pi*0.17)^2 / 10
D_mat <- diag(rep(0.05, 4))
library("ODEnetwork")
lfonet <- ODEnetwork(masses = M_mat, dampers = D_mat, springs = K_mat)
LFOpars <- c("k.1", "k.2", "k.3", "k.4",
             "d.1", "d.2", "d.3", "d.4")
LFObinf <- c(rep(0.2, 4), rep(0.01, 4))
LFObsup <- c(rep(20, 4), rep(0.1, 4))
lfonet <- setState(lfonet, state1 = rep(2, 4), state2 = rep(0, 4))
LFOtimes <- seq(25, 150, by = 2.5)
set.seed(283)
# Warning: The following code might take very long!

LFOres_morris <- ODEmorris(mod = lfonet,
                           pars = LFOpars,
                           times = LFOtimes,
                           binf = LFObinf,
                           bsup = LFObsup,
                           r = 500,
                           design = list(type = "oat", 
                                         levels = 10, grid.jump = 1),
                           scale = TRUE,
                           parallel_eval = TRUE,
                           parallel_eval_ncores = 2)
plot(LFOres_morris, pars_plot = paste0("k.", 1:4), state_plot = "x.2",
     kind = "sep", colors_pars = my_cols)