Particle Markov chain Monte Carlo (PMCMC) algorithm
Usage
pmcmc(
model,
obs_process,
data,
priors,
n_particles,
niter,
theta = NULL,
covmat = NULL,
adaptmix = 0.05,
adaptive = 100,
init_model = NULL,
post_particle = NULL,
record = TRUE,
verbose = getOption("verbose", FALSE)
)
# S4 method for class 'SimInf_model'
pmcmc(
model,
obs_process,
data,
priors,
n_particles,
niter,
theta = NULL,
covmat = NULL,
adaptmix = 0.05,
adaptive = 100,
init_model = NULL,
post_particle = NULL,
record = TRUE,
verbose = getOption("verbose", FALSE)
)
Arguments
- model
The model to simulate data from.
- obs_process
Specification of the stochastic observation process. The
obs_process
can be specified as aformula
if the model contains only one node and there is only one data point for eachtime
indata
. The left hand side of the formula must match a column name in thedata
data.frame and the right hand side of the formula is a character specifying the distribution of the observation process, for example,Iobs ~ poisson(I)
. The following distributions are supported:x ~ binomial(size, prob)
,x ~ poisson(rate)
andx ~ uniform(min, max)
. The observation process can also be a function to evaluate the probability density of the observations given the simulated states. The first argument passed to theobs_process
function is the result from a run of the model and it contains one trajectory with simulated data for a time-point. The second argument to theobs_process
function is adata.frame
containing the rows for the specific time-point that the function is called for. Note that the function must return the log of the density.- data
A
data.frame
holding the time series data.- priors
The priors for the parameters to fit. Each prior is specified with a formula notation, for example,
beta ~ uniform(0, 1)
specifies that beta is uniformly distributed between 0 and 1. Usec()
to provide more than one prior, for example,c(beta ~ uniform(0, 1), gamma ~ normal(10, 1))
. The following distributions are supported:gamma
,lognormal
,normal
anduniform
. All parameters inpriors
must be only in eithergdata
orldata
.- n_particles
An integer with the number of particles (> 1) to use at each timestep.
- niter
An integer specifying the number of iterations to run the PMCMC.
- theta
A named vector of initial values for the parameters of the model. Default is
NULL
, and then these are sampled from the prior distribution(s).- covmat
A named numeric
(npars x npars)
matrix with covariances to use as initial proposal matrix. If left unspecified then defaults todiag((theta/10)^2/npars)
.- adaptmix
Mixing proportion for adaptive proposal. Must be a value between zero and one. Default is
adaptmix = 0.05
.- adaptive
Controls when to start adaptive update. Must be greater or equal to zero. If
adaptive=0
, then adaptive update is not performed. Default isadaptive = 100
.- init_model
FIXME.
- post_particle
An optional function that, if non-NULL, is applied after each completed particle. The function must accept three arguments: 1) an object of
SimInf_pmcmc
with the current state of the fitting process, 2) an objectSimInf_pfilter
with the last particle and one filtered trajectory attached, and 3) an integer with the iteration in the fitting process. This function can be useful to, for example, monitor, save and inspect intermediate results.- record
FIXME
- verbose
prints diagnostic messages when
TRUE
. The default is to retrieve the global optionverbose
and useFALSE
if it is not set. Whenverbose=TRUE
, information is printed every 100 iterations. For pmcmc, it is possible to get information every nth information by specifyingverbose=n
, for example,verbose=1
orverbose=10
.