R/abc_utilities.R

Defines functions ensemble_abc_wrapper create_abc_settings_object graph_Posteriors_All_Parameters

Documented in create_abc_settings_object ensemble_abc_wrapper graph_Posteriors_All_Parameters

#' Wrapper to allow EasyABC functions to run using Ensemble
#'
#' Provides a means of running the ensemble within the EasyABC methods. This
#' method should be stated as the "model" argument of EasyABC methods such as
#' ABC_sequential. Note that as arguments cannot be passed to the model
#' function, these must be created before calling any EasyABC method (see
#' code description for \code{create_abc_settings_object}. The created object
#' most contain the simulation parameters, output measures, an ensemble object
#' to be run, declared as built_ensemble, and whether or not the prior
#' (and generated predictions) should be normalised (normalise). Should this
#' object not exist an error message will be produced.
#' @param x Set of parameter values generated by an EasyABC method
#' @return Ensemble prediction using parameter set in x
#'
#' @export
ensemble_abc_wrapper <- function(x) {

  # Problem with this function is we need to have the parameters, measures,
  # normalise, and ensemble, and we can't pass these in
  # So wrapper object called before this creates an object in the working
  # directory that is then read in here
  if (file.exists(paste(getwd(), "/abc_settings.Rda", sep = ""))) {
    abc_set <- NULL
    load(paste(getwd(), "/abc_settings.Rda", sep = ""))

    params <- as.data.frame(matrix(x, nrow = 1, ncol = length(x)))
    names(params) <- abc_set$parameters

    # Format the parameter input
    prediction <- use_ensemble_to_generate_predictions(
      abc_set$built_ensemble, params, abc_set$parameters, abc_set$measures,
      normalise_values = abc_set$normalise_values,
      normalise_result = abc_set$normalise_result)


    return(prediction)
  } else if(exists("abc_set")) {

    params <- as.data.frame(matrix(x, nrow = 1, ncol = length(x)))
    names(params) <- abc_set$parameters

    # Format the parameter input
    prediction <- use_ensemble_to_generate_predictions(
      abc_set$built_ensemble, params, abc_set$parameters, abc_set$measures,
      normalise_values = abc_set$normalise_values,
      normalise_result = abc_set$normalise_result)

    return(prediction)

  } else {
    message("parameters, measures, and best_ensemble must exist in the workspace,
          declared in an abc_settings object. Run that method first")
  }
}

#' Creates ensemble-specific parameters for ABC analysis
#'
#' The EasyABC model wrapper can only take one parameter input: the parameter
#' values. This is problematic as to generate a prediction for those values,
#' we must provide the names of the simulation parameters and measures, the
#' built ensemble, and whether or not the parameter set and responses have
#' been normalised. To get around that problem, this method creates an object
#' in the working directory that contains these values, and the ensemble abc
#' wrapper provided in spartan can then read these in. Thus, this method
#' MUST be run before using the EasyABC package with the ensemble
#'
#' @param parameters Array containing the names of the parameters for which
#' posterior distributions will be generated
#' @param measures Names of the simulation output responses which the ensemble
#' predicts
#' @param built_ensemble An ensemble object created by spartan
#' @param normalise_values Whether the data provided by the EasyABC algorithm
#' should be normalised (as the ensemble must take data between 0 and 1).
#' More than likley this is TRUE, to ensure the posterior distributions are
#' presented in their correct scale
#' @param normalise_result Whether the results produced in running abc
#' generated parameter sets using the ensemble should be rescaled.
#' @param file_out Whether the settings file should be output to file (TRUE)
#' or as an R object (FALSE)
#' @return Object containing attributes needed for abc analysis
#'
#' @export
create_abc_settings_object <- function(parameters, measures, built_ensemble,
                                       normalise_values = FALSE,
                                       normalise_result = FALSE, file_out = TRUE){

  abc_set <- list("parameters" = parameters, "measures" = measures,
                  "built_ensemble" = built_ensemble,
                  "normalise_values" = normalise_values,
                  "normalise_result" = normalise_result)

  if(file_out)
    save(abc_set, file = paste(getwd(), "/abc_settings.Rda", sep = ""))
  else
    return(abc_set)
}


#' Graph posterior distributions generated for all parameters, to PDF file
#'
#' We provide a means of plotting the produced posterior distribution for all
#' parameters for which the value is being explored. Output to PDF in the
#' working directory
#' @param abc_resultset Result object obtained from the EasyABC package
#' @param parameters Array containing the names of the parameters for which
#' posterior distributions will be generated
#' @param sampleMins Minimum value of the range over which each parameter was
#' explored using ABC
#' @param sampleMaxes Maximum value of the range over which each parameter was
#' explored using ABC
#' @param output_formats File formats in which result graphs should be produced
#'
#' @export
#' @importFrom ggplot2 ggplot aes geom_histogram geom_density scale_x_continuous scale_y_continuous labs theme element_text ggsave geom_point ggtitle geom_abline
#'
graph_Posteriors_All_Parameters <- function(abc_resultset, parameters,
                                            sampleMins, sampleMaxes, output_formats=c("pdf")) {
  for (p in 1:length(parameters)) {

    # Value for each parameter in each column, just need to plot the
    # density of that set
    GRAPHTITLE <- "Parameter Value Search Using \n Approximate Bayesian
      Computation"

    for(out_type in output_formats)
    {
      GRAPH <- ggplot(data.frame(abc_resultset$param[, p]),
                      aes(abc_resultset$param[, p])) +
        geom_density() +
        scale_x_continuous(limits = c(sampleMins[p], sampleMaxes[p])) +
        labs(x = paste(parameters[p], " Value", sep = ""),
             y = "Density", title = GRAPHTITLE,
             subtitle = paste("Parameter: ", parameters[p], sep = "")) +
        theme(axis.title = element_text(size = 11),
              axis.text = element_text(size = 11),
              plot.title = element_text(size = 11, hjust = 0.5),
              plot.subtitle = element_text(size = 11, hjust = 0.5))

      ggsave(paste0(parameters[p], ".",out_type),
             plot = GRAPH, width = 4, height = 4)
    }
  }
}

Try the spartan package in your browser

Any scripts or data that you put into this service are public.

spartan documentation built on May 2, 2019, 9:39 a.m.