Main wrapper function to start MCMCs, particle MCMCs and SMCs
either a BayesianSetup (see
sampling algorithm to be run. Default is DEzs. Options are "Metropolis", "AM", "DR", "DRAM", "DE", "DEzs", "DREAM", "DREAMzs", "SMC". For details see the help of the individual functions.
list with settings for each sampler (see help of sampler for details). If a setting is not provided, defaults (see
The runMCMC function can be started with either one of a) an object of class BayesianSetup with prior and likelihood function (recommended, see
createBayesianSetup), b) a log posterior or other target function, or c) an object of class BayesianOutput created by runMCMC. The latter allows to continue a previous MCMC run. If a bayesianSetup is provided, check if appropriate parallization options are used - many samplers can make use of parallelization if this option is activated when the class is created.
For details about the different MCMC samplers, make sure you have read the Vignette (run vignette("BayesianTools", package="BayesianTools"). Also, see
Metropolis for Metropolis based samplers,
DEzs for standard differential evolution samplers,
DREAMzs for DREAM sampler,
Twalk for the Twalk sampler, and
smcSampler for rejection and Sequential Monte Carlo sampling.
The samplers "AM", "DR", and "DRAM" are special cases of the "Metropolis" sampler and are shortcuts for predefined settings ("AM": adapt=TRUE; "DR": DRlevels=2; "DRAM": adapt=True, DRlevels=2).
The settings list allows to change the settings for the MCMC samplers and some other options. For the MCMC sampler settings, see their help files. Global options that apply for all MCMC samplers are: iterations (number of MCMC iterations), and nrChains (number of chains to run). Note that running several chains is not done in parallel, so if time is an issue it will be better to run the MCMCs individually and then combine them via
createMcmcSamplerList into one joint object.
Startvalues: all samplers allow to provide explicit startvalues. If startvalues are not provided, they are sampled from the prior. Usually, this is a good choice, so don't feel compelled to provide startvalues.
Note that DE and DREAM variants as well as SMC and T-walk require a population to start, which should be provided as a matrix. Default (NULL) sets the population size for DE to 3 x dimensions of parameters, for DREAM to 2 x dimensions of parameters and for DEzs and DREAMzs to three, sampled from the prior. Note also that the zs variants of DE and DREAM require two populations, the current population and the z matrix (a kind of memory) - if you want to set both, provide a list with startvalue$X and startvalue$Z.
setting startValue for sampling with nrChains > 1 : if you want to provide different start values for the different chains, provide them as a list
The function returns an object of class mcmcSampler (if one chain is run) or mcmcSamplerList. Both have the superclass bayesianOutput. It is possible to extract the samples as a coda object or matrix with
It is also possible to summarize the posterior as a new prior via
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
## Generate a test likelihood function. ll <- generateTestDensityMultiNormal(sigma = "no correlation") ## Create a BayesianSetup object from the likelihood ## is the recommended way of using the runMCMC() function. bayesianSetup <- createBayesianSetup(likelihood = ll, lower = rep(-10, 3), upper = rep(10, 3)) ## Finally we can run the sampler and have a look settings = list(iterations = 1000, adapt = FALSE) out <- runMCMC(bayesianSetup = bayesianSetup, sampler = "Metropolis", settings = settings) ## out is of class bayesianOutput. There are various standard functions # implemented for this output plot(out) correlationPlot(out) marginalPlot(out) summary(out) ## additionally, you can return the sample as a coda object, and make use of the coda functions # for plotting and analysis codaObject = getSample(out, start = 500, coda = TRUE)
sh: 1: cannot create /dev/null: Permission denied sh: 1: cannot create /dev/null: Permission denied BT runMCMC: trying to find optimal start and covariance values BT runMCMC: Optimization finished, setting startValues to 0.000475896965443028 -0.000587200694203585 0.000248667831606472 - Setting covariance to 0.999943540247535 2.94245811258299e-05 7.51586444568408e-05 2.94245811258299e-05 0.999970969717136 4.68662602402109e-06 7.51586444568408e-05 4.68662602402109e-06 0.999756376304926 Running Metropolis-MCMC, chain 1 iteration 100 of 1000 . Current logp: -12.92031 Please wait! Running Metropolis-MCMC, chain 1 iteration 200 of 1000 . Current logp: -14.15851 Please wait! Running Metropolis-MCMC, chain 1 iteration 300 of 1000 . Current logp: -12.82188 Please wait! Running Metropolis-MCMC, chain 1 iteration 400 of 1000 . Current logp: -12.28031 Please wait! Running Metropolis-MCMC, chain 1 iteration 500 of 1000 . Current logp: -14.38009 Please wait! Running Metropolis-MCMC, chain 1 iteration 600 of 1000 . Current logp: -12.45922 Please wait! Running Metropolis-MCMC, chain 1 iteration 700 of 1000 . Current logp: -14.08274 Please wait! Running Metropolis-MCMC, chain 1 iteration 800 of 1000 . Current logp: -12.05949 Please wait! Running Metropolis-MCMC, chain 1 iteration 900 of 1000 . Current logp: -15.0019 Please wait! Running Metropolis-MCMC, chain 1 iteration 1000 of 1000 . Current logp: -14.59789 Please wait! runMCMC terminated after 1.469seconds Package 'sm', version 2.2-5.4: type help(sm) for summary information # # # # # # # # # # # # # # # # # # # # # # # # # ## MCMC chain summary ## # # # # # # # # # # # # # # # # # # # # # # # # # # MCMC sampler: Metropolis # Nr. Chains: 1 # Iterations per chain: 1000 # Rejection rate: 0.652 # Effective sample size: 129 # Runtime: 1.469 sec. # Parameter MAP 2.5% median 97.5% # par 1 : 0.000 -2.170 0.010 1.792 # par 2 : -0.001 -1.787 -0.038 2.064 # par 3 : 0.000 -1.623 -0.019 1.864 ## DIC: 10.904 ## Convergence Gelman Rubin multivariate psrf: Only one chain; convergence cannot be determined! ## Correlations par 1 par 2 par 3 par 1 1.000 -0.081 -0.037 par 2 -0.081 1.000 0.044 par 3 -0.037 0.044 1.000
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.