SABC | R Documentation |

Algorithms for the Simulated Annealing approach to Approximate Bayesian Computation (SABC).

SABC(r.model, r.prior, d.prior, n.sample, eps.init, iter.max, v=ifelse(method=="informative",0.4,1.2), beta=0.8, delta=0.1, resample=5*n.sample, verbose=n.sample, method="noninformative", adaptjump=TRUE, summarystats=FALSE, y=NULL, f.summarystats=NULL, ...)

`r.model` |
Function that returns either a random sample from the likelihood or a scalar distance between such a sample and the data. The first argument must be the parameter vector. |

`r.prior` |
Function that returns a random sample from the prior. |

`d.prior` |
Function that returns the density of the prior distribution. |

`n.sample` |
Size of the ensemble. |

`eps.init` |
Initial tolerance or temperature. |

`iter.max` |
Total number of simulations from the likelihood. |

`v` |
Tuning parameter that governs the annealing speed. Defaults to 1.2, for the |

`beta` |
Tuning parameter that governs the mixing in parameter space. Defaults to 0.8. |

`delta` |
Tuning parameter for the resampling steps. Defaults to 0.1. |

`resample` |
Number of accepted particle updates after which a resampling step is performed. Defaults to 5* |

`verbose` |
Shows the iteration progress each |

`adaptjump` |
Whether to adapt covariance of jump distribution. Default is TRUE. |

`method` |
Argument to select algorithm. Accepts |

`summarystats` |
Whether summary statistics shall be calculated (semi-) automatically. Defaults to FALSE. |

`y` |
Data vector. Needs to be provided if either |

`f.summarystats` |
If |

`...` |
further arguments passed to |

SABC defines a class of algorithms for particle ABC that are inspired by Simulated Annealing. Unlike other algorithms, this class is not based on importance sampling, and hence does not suffer from a loss of effective sample size due to re-sampling. The approach is presented in detail in Albert, Kuensch, and Scheidegger (2014; see references).

This package implements two versions of SABC algorithms, for the cases of a non-informative or an informative prior. These are described in detail in the paper. The algorithms can be selected using the `method`

argument: `method=noninformative`

or `method=informative`

.
In the informative case, the algorithm corrects for the bias caused by an over- or under-representation of the prior.

The argument `adaptjump`

allows a choice of whether to adapt the covariance of the jump distribution. Default is TRUE.

Furthermore, the package allows for three different ways of using the data.
If `y`

is not provided, the algorithm expects `r.model`

to return a scalar measuring the distance between a random sample from the likelihood and the data.
If `y`

is provided and `summarystats = FALSE`

, the algorithm expects `r.model`

to return a random sample from the likelihood and uses the relative sum of squares to measure the distances between `y`

and random likelihood samples.
If `summarystats = TRUE`

the algorithm calculates summary statistics semi-automatically, as described in detail in the paper by Fearnhead et al. (2012; see references).
The summary statistics are calculated by means of a linear regression applied to a sample from the prior and the image of `f.summarystats`

of an associated sample from the likelihood.

Returns a list with the following components:

`E` |
Matrix with ensemble of samples. |

`P` |
Matrix with prior ensemble of samples. |

`eps` |
Value of tolerance (temperature) at final iteration. |

`ESS` |
Effective sample size, due to final bias correction ( |

Carlo Albert <carlo.albert@eawag.ch>, Andreas Scheidegger, Tobia Fasciati. Package initially compiled by Lukas M. Weber.

C. Albert, H. R. Kuensch and A. Scheidegger, Statistics and Computing 0960-3174 (2014), arXiv:1208.2157, *A Simulated Annealing Approach to Approximate Bayes Computations*.

P. Fearnhead and D. Prangle, J. R. Statist. Soc. B 74 (2012), *Constructing summary statistics for approximate Bayesian computation: semi-automatic approximate Bayesian computation*.

## Not run: ## Example for "noninformative" case # Prior is uniform on [-10,10] d.prior <- function(par) dunif(par,-10,10) r.prior <- function() runif(1,-10,10) # Model is the sum of two normal distributions. Return distance to observation 0: f.dist <- function(par) return( abs(rnorm( 1 , par , ifelse(runif(1)<0.5,1,0.1 ) ))) # Run algorithm ("noninformative" case) res <- SABC(f.dist,r.prior,d.prior,n.sample=500,eps.init=2,iter.max=50000) ## End(Not run) ## Not run: # Histogram of results hist(res$E[,1],breaks=200) ## End(Not run)

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.