ABC_emulation: Rejection sampling scheme for ABC using an emulator

View source: R/ABC_emulation.R

ABC_emulationR Documentation

Rejection sampling scheme for ABC using an emulator

Description

This function launches a series of nb_design_pts model simulations with model parameters drawn in the prior distribution specified in prior_matrix, build an emulator with these computed design points and then launches a series of nb_simul emulator simulations.

Usage

  ABC_emulation(model, prior,  nb_design_pts, nb_simul, prior_test=NULL,
  summary_stat_target=NULL, emulator_span = 50, tol=NULL, use_seed=FALSE,
  seed_count=0, n_cluster=1, verbose=FALSE, progress_bar=FALSE)

Arguments

model

a R function implementing the model to be simulated. It must take as arguments a vector of model parameter values and it must return a vector of summary statistics. When using the option use_seed=TRUE, model must take as arguments a vector containing a seed value and the model parameter values. A tutorial is provided in the package's vignette to dynamically link a binary code to a R function. Users may alternatively wish to wrap their binary executables using the provided functions binary_model and binary_model_cluster. The use of these functions is associated with slightly different constraints on the design of the binary code (see binary_model and binary_model_cluster).

prior

a list of prior information. Each element of the list corresponds to a model parameter. The list element must be a vector whose first argument determines the type of prior distribution: possible values are "unif" for a uniform distribution on a segment, "normal" for a normal distribution, "lognormal" for a lognormal distribution or "exponential" for an exponential distribution. The following arguments of the list elements contain the characteritiscs of the prior distribution chosen: for "unif", two numbers must be given: the minimum and maximum values of the uniform distribution; for "normal", two numbers must be given: the mean and standard deviation of the normal distribution; for "lognormal", two numbers must be given: the mean and standard deviation on the log scale of the lognormal distribution; for "exponential", one number must be given: the rate of the exponential distribution. User-defined prior distributions can also be provided. See the vignette for additional information on this topic.

nb_design_pts

a positive integer equal to the desired number of simulations of the model used to build the emulator.

nb_simul

a positive integer equal to the desired number of simulations of the emulator.

prior_test

a string expressing the constraints between model parameters. This expression will be evaluated as a logical expression, you can use all the logical operators including "<", ">", ... Each parameter should be designated with "X1", "X2", ... in the same order as in the prior definition. If not provided, no constraint will be applied.

summary_stat_target

a vector containing the targeted (observed) summary statistics. If not provided, ABC_rejection only launches the simulations and outputs the simulation results.

emulator_span

a positive number, the number of design points selected for the local regression. 50 by default.

tol

tolerance, a strictly positive number (between 0 and 1) indicating the proportion of simulations retained nearest the targeted summary statistics.

use_seed

logical. If FALSE (default), ABC_rejection provides as input to the function model a vector containing the model parameters used for the simulation. If TRUE, ABC_rejection provides as input to the function model a vector containing an integer seed value and the model parameters used for the simulation. In this last case, the seed value should be used by model to initialize its pseudo-random number generators (if model is stochastic).

seed_count

a positive integer, the initial seed value provided to the function model (if use_seed=TRUE). This value is incremented by 1 at each call of the function model.

n_cluster

a positive integer. If larger than 1 (the default value), ABC_rejection will launch model simulations in parallel on n_cluster cores of the computer.

verbose

logical. FALSE by default. If TRUE, ABC_rejection writes in the current directory intermediary results at the end of each step of the algorithm in the file "output". These outputs have a matrix format, in wich each raw is a different simulation, the first columns are the parameters used for this simulation, and the last columns are the summary statistics of this simulation.

progress_bar

logical, FALSE by default. If TRUE, ABC_rejection will output a bar of progression with the estimated remaining computing time. Option not available with multiple cores.

Value

The returned value is a list containing the following components:

param

The model parameters used in the model simulations.

stats

The summary statistics obtained at the end of the model simulations.

weights

The weights of the different model simulations. In the standard rejection scheme, all model simulations have the same weights.

stats_normalization

The standard deviation of the summary statistics across the model simulations.

nsim

The number of model simulations performed.

nrec

The number of retained simulations (if targeted summary statistics are provided).

computime

The computing time to perform the simulations.

Author(s)

Franck Jabot, Thierry Faure and Nicolas Dumoulin

References

Jabot, F., Lagarrigues G., Courbaud B., Dumoulin N. (2015). A comparison of emulation methods for Approximate Bayesian Computation. To be published.

See Also

binary_model, binary_model_cluster, ABC_sequential, ABC_mcmc

Examples

 ## Not run:  

    ##### EXAMPLE 1 #####
    #####################

    ## the model is a C++ function packed into a R function -- the option 'use_seed'
    ##  must be turned to TRUE.
    trait_prior=list(c("unif",3,5),c("unif",-2.3,1.6),c("unif",-25,125),c("unif",-0.7,3.2))
    trait_prior

    ## only launching simulations with parameters drawn in the prior distributions
    ABC_emul = ABC_emulation(model=trait_model, prior=trait_prior,
      nb_design_pts=10, nb_simul=300, use_seed=TRUE, progress=TRUE)
    ABC_emul

    ## launching simulations with parameters drawn in the prior distributions and performing
    # the rejection step
    sum_stat_obs=c(100,2.5,20,30000)
    ABC_emul = ABC_emulation(model=trait_model, prior=trait_prior, tol=0.2, nb_design_pts=10,
      nb_simul=100, summary_stat_target=sum_stat_obs, use_seed=TRUE, progress=TRUE)
    ABC_emul

 
## End(Not run)

EasyABC documentation built on Jan. 6, 2023, 1:14 a.m.