Construct parameters for use with
pmcmc()
. This creates a utility object that is used
internally to work with parameters. Most users only need to
construct this object, but see the examples for how it can be
used.
Unless you have a very simple model, it is highly unlikely that the parameters that you are interested in performing inference on are the same as the parameters that you might need to initialise your model.
Due to the nature of mcmc and other inference algorithms, the
general assumption is that the inference parameters will be a
simple vector of real values; here each of the parameters
elements corresponds to one of these. The proposal matrix maps
one vector to another via a simple multivariate-gaussian kernel.
On the other hand, dust models can take a named list of arbitrary data as their input parameters (see dust::dust_generator). These might include:
things that are not parameters at all from the perspective of the inference - for example some quantity that you might vary depending on the region/species/etc you're running the model for but that you are not fitting.
non-scalar quantities that are directly derived from some parameters that you are fitting. As an example of this, in sircovid, a transmission model of COVID, we take a number of "contact rates" which apply at different points in time, and generate from this an interpolated series of contact rates per time step (a very long vector). Other users have needed to generate equilibrium solutions to parts of their model and used these at initialisation.
arbitrary complex inputs to the model, for example weather data, demographic matrices, population contact rate matrices etc. These are all "parameters" from the perspective of a dust model but not at all from the perspective of the inference process.
To allow for this in a flexible way, mcstate allows a "transform"
function, the transform
argument to the constructor. This
function maps a named numeric vector of inference parameters to
whatever you need for your dust model. The default value for
this function is as.list which just converts the named vector
to a named list, which works well in the example cases here.
When providing a transformation function, you may want to provide
a "closure" rather than a top-level function. This way you can
bind additional data into your function. For example, suppose
that you want to use some demographic matrix m
in your model,
and perform inference on parameters a
and b
you might write
and pass this into mcstate::pmcmc_parameters$new
, providing
parameter definitions only for a
and b
. See the examples
for full working of this.
new()
Create the pmcmc_parameters object
pmcmc_parameters$new(parameters, proposal, transform = NULL)
parameters
A list
of
pmcmc_parameter objects, each of which describe a
single parameter in your model. If parameters
is named, then
these names must match the $name
element of each parameter is
used (this is verified).
proposal
A square proposal distribution corresponding to the
variance-covariance matrix of a multivariate gaussian distribution
used to generate new parameters. It must have the same number of
rows and columns as there are elements in parameters
, and if
named the names must correspond exactly to the names in
parameters
. Because it corresponds to a variance-covariance
matrix it must be symmetric and positive definite.
transform
An optional transformation function to apply to your parameter vector immediately before passing it to the model function. If not given, then as.list is used, as dust models require this. However, if t you need to generate derived parameters from those being actively sampled you can do arbitrary transformations here.
mean()
Return the estimate of the mean of the parameters, as set when created (this is not updated by any fitting!)
summary()
Return a data.frame
with information about
parameters (name, min, max, and integer).
propose()
Propose a new parameter vector given a current parameter
vector. This proposes a new parameter vector given your current
vector and the variance-covariance matrix of your proposal
kernel, rounds any integer values, and reflects bounded
parameters until they lie within min
:max
.
theta
a parameter vector in the same order as your
parameters were defined in (see $names()
for that order.
scale
an optional scaling factor to apply to the proposal distribution. This may be useful in sampling starting points. The parameter is equivalent to a multiplicative factor applied to the variance covariance matrix.
vcv
A variance covariance matrix of the correct size, overriding the proposal matrix built into the parameters object. This will be slightly less efficient but allow a different proposal matrix to be used (e.g., during an adaptive MCMC)
fix()
Set some parameters to fixed values. Use this to reduce the dimensionality of your system.
# Construct an object with two parameters:
pars <- mcstate::pmcmc_parameters$new(
list(mcstate::pmcmc_parameter("a", 0.1, min = 0, max = 1,
prior = function(a) log(a)),
mcstate::pmcmc_parameter("b", 0, prior = dnorm)),
matrix(c(1, 0.5, 0.5, 2), 2, 2))
# Initial parameters
p <- pars$initial()
p
#> a b
#> 0.1 0.0
# Propose a new parameter point
pars$propose(p)
#> a b
#> 0.8378669 0.1972515
# Information about parameters:
pars$names()
#> [1] "a" "b"
pars$summary()
#> name min max discrete integer
#> 1 a 0 1 FALSE FALSE
#> 2 b -Inf Inf FALSE FALSE
# Compute prior
pars$prior(p)
#> [1] -1.903643
# Transform data for your model
pars$model(p)
#> $a
#> [1] 0.1
#>
#> $b
#> [1] 0
#>
# Above we describe a nontrivial transformation function using a closure
make_transform <- function(m) {
function(theta) {
c(list(m = m), as.list(theta))
}
}
# Suppose this is our demographic matrix (note here that the name
# need not match that used in the transform)
demographic_matrix <- diag(4)
# Construct the parameters as above, but this time passing in the
# function that make_transform returns
pars <- mcstate::pmcmc_parameters$new(
list(mcstate::pmcmc_parameter("a", 0.1, min = 0, max = 1,
prior = function(a) log(a)),
mcstate::pmcmc_parameter("b", 0, prior = dnorm)),
matrix(c(1, 0.5, 0.5, 2), 2, 2),
make_transform(demographic_matrix))
# Now, as above we start from a position in terms of a and b only:
pars$initial()
#> a b
#> 0.1 0.0
# But when prepared for the model, our matrix will be set up
pars$model(pars$initial())
#> $m
#> [,1] [,2] [,3] [,4]
#> [1,] 1 0 0 0
#> [2,] 0 1 0 0
#> [3,] 0 0 1 0
#> [4,] 0 0 0 1
#>
#> $a
#> [1] 0.1
#>
#> $b
#> [1] 0
#>