R/e2e_plot_opt_diagnostics.R

In StrathE2E2: End-to-End Marine Food Web Model

#
# e2e_plot_opt_diagnostics.R
#
#' Plot diagnostic data on the performance of parameter optimization functions.
#'
#' Creates diagnostic plots from outputs of the functions for optimizing ecology, harvest ratio, and fishing activity parameters: e2e_optimize_eco(), e2e_optimize_hr(), and e2e_optimize_act() respectively. 
#' 
#' The function takes the history of parameter values proposed during each simulated annealing process, expresses these relative to the initial value of each parameter, generates qunatiles (0.005, 0.25, 0.5, 0.75 and 0.995) of the distribution for each parameter, and plots the results as box-and-whisker diagrams.
#' The relativity to initial values is expressed as (proposed - initial)/initial so the initial is represented by a relative value of zero. For each parameter, the final accepted value is over-plotted onto the box-and-whisker as a red bar. The resulting diagram, with a separate box-and-whisker for each parameter
#' shows the extent to which each parameter has migrated from its initial value en-route to the final accepted state.
#'
#' The format of the parameter diagram is the same for each type of optimization scheme: ecology, harvest ratio, and fishing activity parameters. However, in the case where activity parameters are
#' optimized to zonal harvest ratios rather than ecosystem state data, the diagram has an additional panel showing the distribution of proposal harvest ratios relative to the target, and (in red) the harvest 
#' ratios corresponding to the final accepted activity parmeters (relative to the targets).
#'
#' Arguments for this function permit the input data to be drawn from an existing data object generated by the various optimization functions, previously genearted csv files, or example data
#' provided with the package for versions of the internal North Sea models.
#'
#' Documentation in a dataframe format on each of the classes of parameters in the model can be obtained with the function e2e_get_parmdoc(). This provides a key to the abbreviated parameter names especially in the diagnostic plots for e2e_optimize_eco().
#'
#' As well as drawing the plot the function returns a list object containing a) an array of the quantiles of the proposal distributikn for each parameter, and b) an array of one row, of the final accepted parameter set from the procedure.
#'
#' @param model R-list object defining the baseline model configuration used to generate the data and compiled by the e2e_read() function.
#' @param selection Text string from a list identifying which type of optimization procedure generated the data to be plotted. Select from: "ECO", "HR", "ACT", corresponding to the functions e2e_optimize_eco(), e2e_optimize_hr(), and e2e_optimize_act(). Remember to include the phrase within "" quotes.
#' @param fitted.to Specific to the case where selection="ACT"; text string from a list identifying which version of activity optimization procedure generated the data to be plotted. Select from: "ECO", "HR", corresponding to maximising the likelihod of ecosystem state data, or zonal harvest ratios respectively. Default="". Remember to include the phrase within "" quotes.
#' @param use.saved Logical. If TRUE use data from a prior user-defined run held as csv files data in the current results folder (default=FALSE).
#' @param use.example Logical. If TRUE use pre-computed example data from the internal North Sea model rather than user-generated data (default=FALSE).
#' @param results R-list object of output from an optimization process generated by the e2e_optimize_eco(), e2e_optimize_hr(), or e2e_optimize_act() function. Only needed if use.saved1=FALSE and use.example1=FALSE. (Default=NULL).
#'
#' @return List object of results from which the plot is created, graphical display in a new graphics window.
#'
#' @seealso \code{\link{e2e_read}}, \code{\link{e2e_get_parmdoc}}, \code{\link{e2e_optimize_eco}}, \code{\link{e2e_optimize_hr}}, \code{\link{e2e_optimize_act}}
#'
#' @importFrom graphics text
#'
#' @export
#'
#' @examples
#' # The examples provided here are illustration of how to generate diagnostics plots
#' # from optimization runs of the model. Optimization runs are very time consuming so
#' # the examples only involve a bare minimum of model runs and are not realistic.
#' # Alternatively, data in the package StrathE2E2examples can be used to generate 
#' # example plots.
#'
#' # --------------------------------------------------------------------------
#'
#' \donttest{
#' # Load the 1970-1999 version of the North Sea model supplied with the package and generate
#' # a quick test data object with only 8 itereations and running the model for 3 years.
#' # More realistic would be at least 500 iterations and running for 50 years.
#' # Even so this example will take a few minutes to run:
#'     model<-e2e_read(model.name="North_Sea",
#'                     model.variant="1970-1999",
#'                     model.ident="test")
#' # This model is already optimized to the observed ecosystem data supplied with the package
#' # Perturb the temperature driving to knock the model away from its maximum likelihood
#' # state relative to the target data to illustrate the performance:
#' # add 3 degC to upper layer offshore temperatures:
#'     model$data$physics.drivers$so_temp <- model$data$physics.drivers$so_temp+3
#' # add 3 degC to inshore temperatures:
#'     model$data$physics.drivers$si_temp <- model$data$physics.drivers$si_temp+3
#' # add 3 degC to lower layer offshore temperatures:
#'     model$data$physics.drivers$d_temp  <- model$data$physics.drivers$d_temp+3
#'     test_run  <- e2e_optimize_eco(model, nyears=3, n_iter=8, start_temperature=0.4,
#'                                   csv.output=FALSE)
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="ECO",results=test_run)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#'     str(plot_data,max.level=1)  # show the structure of the list object plot_data
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' # Or... plot example date supplied with the package showing some data generated during
#' # the process of optimizing the North Sea model:
#' # This example requires the Strathe2E2examples supplementary data package.
#' if(require(StrathE2E2examples)){
#'     model     <- e2e_read(model.name="North_Sea", model.variant="1970-1999")
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="ECO",use.example=TRUE)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#' # Example data are only available for the 1970-1999 variant of the North Sea model
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' \donttest{
#' # Same for harvest ratio optimization...
#'     model<-e2e_read(model.name="North_Sea",
#'                     model.variant="2003-2013",
#'                     model.ident="test")
#' # This model is already optimized to the observed ecosystem data supplied with the package
#' # Perturb the temperature driving to knock the model away from its maximum likelihood
#' # state relative to the target data to illustrate the performance o fthe process:
#' # add 3 degC to upper layer offshore temperatures:
#'     model$data$physics.drivers$so_temp <- model$data$physics.drivers$so_temp+3
#' # add 3 degC to inshore temperatures:
#'     model$data$physics.drivers$si_temp <- model$data$physics.drivers$si_temp+3
#' # add 3 degC to lower layer offshore temperatures:
#'     model$data$physics.drivers$d_temp  <- model$data$physics.drivers$d_temp+3
#'     test_run  <- e2e_optimize_hr(model, nyears=3, n_iter=8, start_temperature=0.4,
#'                                  csv.output=FALSE)
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="HR",results=test_run)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#'     str(plot_data,max.level=1)  # show the structure of the list object plot_data
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' # Or... plot example date supplied with the package showing some data generated
#' # during the process of optimizing the North Sea model:
#' # This example requires the Strathe2E2examples supplementary data package.
#' if(require(StrathE2E2examples)){
#'     model     <- e2e_read(model.name="North_Sea", model.variant="2003-2013")
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="HR",use.example=TRUE)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#' # Example data are only available for the 2003-2013 variant of the North Sea model
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' \donttest{
#' # For activity rate optimization relative to ecosystem data:
#'     model<-e2e_read(model.name="North_Sea",
#'                     model.variant="1970-1999",
#'                     model.ident="test")
#' # This model is already optimized to the observed ecosystem data supplied with the package,
#' # but not by optimizing gear activity rates
#' # The e2e_optimize_eco() function was used in this case.
#' # Perturb the temperature driving to knock the model away from its maximum likelihood
#' # state relative to the target data:
#' # add 3 degC to upper layer offshore temperatures:
#'     model$data$physics.drivers$so_temp <- model$data$physics.drivers$so_temp+3
#' # add 3 degC to inshore temperatures:
#'     model$data$physics.drivers$si_temp <- model$data$physics.drivers$si_temp+3
#' # add 3 degC to lower layer offshore temperatures:
#'     model$data$physics.drivers$d_temp  <- model$data$physics.drivers$d_temp+3
#'     test_run  <- e2e_optimize_act(model, selection="ECO", n_iter=8, start_temperature=0.4,
#'                                   cooling=0.975, csv.output=FALSE, nyears=3)
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="ACT",
#'                                           fitted.to="ECO",results=test_run)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#'     str(plot_data,max.level=1)  # show the structure of the list object plot_data
#'
#' # There are no example data available in the package for this function
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' \donttest{
#' # For activity rate optimization relative to zonal harvest ratios:
#'     model<-e2e_read(model.name="North_Sea",
#'                     model.variant="1970-1999",
#'                     model.ident="test")
#' # Activity rates in this model are already optimized to the target harvest ratios supplied with
#' # the package but we woud not expect to recover these values in this short demonstration run
#'     test_run  <- e2e_optimize_act(model, selection="HR", n_iter=30, start_temperature=1.0,
#'                                   cooling=0.985, csv.output=FALSE, n_traj=5 )
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="ACT",fitted.to="HR",
#'                                           results=test_run)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#'     str(plot_data,max.level=1)  # show the structure of the list object plot_data
#' }
#'
#' # --------------------------------------------------------------------------
#'
#' # Or... plot example date supplied with the package showing some data generated during
#' # the process of optimizing the North Sea model:
#' # This example requires the Strathe2E2examples supplementary data package.
#' if(require(StrathE2E2examples)){
#'     model     <- e2e_read(model.name="North_Sea", model.variant="1970-1999")
#'     plot_data <- e2e_plot_opt_diagnostics(model,selection="ACT",fitted.to="HR",
#'                                           use.example=TRUE)
#' # Red bars shows the final accepted values of each parameters, boxes and whiskers
#' # show the dirtribution of paramater values explored, relative to the starting values
#' # Example data are only available for the 1970-1999 variant of the North Sea model
#' }
#'
#' # --------------------------------------------------------------------------
#'
#
# ---------------------------------------------------------------------
# |                                                                   |
# | Author: Mike Heath                                                |
# | Department of Mathematics and Statistics                          |
# | University of Strathclyde, Glasgow                                |
# |                                                                   |
# | Date of this version: May 2020                                    |
# |                                                                   |
# ---------------------------------------------------------------------

e2e_plot_opt_diagnostics <- function(model,selection="",fitted.to="",use.saved=FALSE, use.example=FALSE, results=NULL) {

   oo <- options()
   on.exit(options(oo))

	if (use.example == TRUE) hasExamples()

output_data<-NULL

allowed.list_selection<-c("ECO","HR","ACT")
if(selection %in% allowed.list_selection==FALSE){
	message("The argument 'selection' must be one of the following...")
	message("  'ECO' = results from optimizing ecology parameters")
	message("  'HR'  = results from optimizing harvest ratio multipliers")
	message("  'ACT' = results from optimizing gear activity multipliers")
	stop(" unknown selection '", selection, "' !\n")
}

allowed.list_fitted.to<-c("ECO","HR")
if(selection=="ACT" & fitted.to %in% allowed.list_fitted.to==FALSE){
	message("With 'selection' = 'ACT', the argument 'fitted.to' must be one of the following...")
	message("  'ECO' = results from optimizing to observed ecosystem data")
	message("  'HR'  = results from optimizing to known harvest ratios")
	stop(" unknown setting '", fitted.to, "' !\n")
}


if(selection=="ECO"){
output_data<-plot_eco_hr_act_eco_optimization_diagnostics(model=model, selection=selection, use.saved=use.saved, use.example=use.example, results=results)

} else if(selection=="HR"){
output_data<-plot_eco_hr_act_eco_optimization_diagnostics(model=model, selection=selection, use.saved=use.saved, use.example=use.example, results=results)

} else if(selection=="ACT" & fitted.to=="ECO"){
output_data<-plot_eco_hr_act_eco_optimization_diagnostics(model=model, selection=selection, use.saved=use.saved, use.example=use.example, results=results)

} else if(selection=="ACT" & fitted.to=="HR"){
output_data<-plot_act_hr_optimization_diagnostics(model=model, use.saved=use.saved, use.example=use.example, results=results)

} else {
	stop("Error: unknown 'selection' and or 'fitted.to' values '", selection, "' and '", fitted.to, "' !\n")
}


return(output_data)

}

#-----------------------------------------------------------------