Bayesian model comparison using the Hyvärinen score.
[more details in S. Shao, P.E. Jacob, J. Ding, and V. Tarokh (2017), available here (or alternatively there)]
This package provides functions that compute the Hyvärinen score (and the log-evidence as an aside) of a Bayesian model. This is achieved by using either SMC [e.g. Chopin (2002) and Del Moral, Doucet, Jasra (2006)] or SMC2 [cf. Chopin, Jacob, Papaspiliopoulos (2013)] depending on whether the likelihood can be evaluated.
- Computing the Hyvärinen score
- Defining a model
- Setting algorithmic parameters
- Output description
- Particle filter
- Examples
After installing and loading the package, computing the Hyvärinen score (abbrv. H-score) of a model is done by calling
hscore(observations, model, algorithmic_parameters)
where the inputs consist of
observations
: sequence of dimY-dimensional observations Y1 , ... , YT encoded as a dimY by Tmatrix
model
: Bayesian model encoded as alist
(see Defining a model)algorithmic_parameters
: algorithmic parameters provided as alist
(see Setting algorithmic parameters)
The function hscore
is a wrapper that either calls smc(observations, model, algorithmic_parameters)
or smc2(observations, model, algorithmic_parameters)
depending on whether the model
specifies the likelihood or not.
The output is detailed below (see Output description).
The model
needs to be provided as a list
. A description of the required fields is provided in inst/_template_model.R
. Missing fields may be automatically filled-in via the function set_default_model
defined in R/util_default.R
.
The algorithmic_parameters
need to be provided as a list
. A complete description of the required fields along with their default values can be found in the function set_default_algorithmic_parameters
defined in R/util_default.R
.
The output of hscore
, smc
, or smc2
is a list
. Depending on the specified algorithmic_parameters
, some of its fields may be set to NULL
. In its most exhaustive form, the output will contain the following objects:
thetas
: last set of particles thetas (dimtheta by Nthetamatrix
)normw
: corresponding normalized weights (vector
of length Ntheta)byproducts
orPFs
: correspondinglist
of byproducts (e.g. particle filters forsmc2
, see Particle filter)logtargetdensities
: corresponding evaluations of target log-densities (vector
of length Ntheta)thetas_history
:list
of successive sets of particles thetas (onematrix
per timestep, starting with the prior, so the length of thelist
is T+1)normw_history
:list
of successive normalized weights (i.e.list
ofvector
)logtargetdensities_history
:list
of successive target log-densities evaluations (i.e.list
ofvector
)byproducts_history
ofPF_history
:list
of successive byproducts or particle filters (i.e.list
oflist
)logevidence
: cumulative logevidence (vector
of length T)hscore
: cumulative H-score using Fisher/Louis type identities (vector
of length T)hscoreDDE
: cumulative H-score using kernel density estimation (vector
of length T)ESS
: successive ESS (vector
of length T)rejuvenation_times
: successive times when rejuvenation occured (vector
of random length less than T)rejuvenation_rate
: associated acceptance rates (vector
of random length less than T)method
: method called (string
equal to either"SMC"
or"SMC2"
)algorithmic_parameters
:list
of algorithmic parameters used
The particle filter is implemented in R/conditional_particle_filter.R
. The function conditional_particle_filter
uses the bootstrap particle filter. It takes an optional argument path
when conditioning on a particular path is needed. The output is a list
containing the number of x-particles (Nx
), the last set of x-particles (X
), their respective normalized weights (xnormW
), an estimator of the log-likelihood (log_p_y_hat
), the last incremental log-likelihood (incremental_ll
), a tree encoding the paths and ancestral lineages (tree
), and a realized path sampled from the final sets of paths (path
).
Some examples are provided in the folder inst
. These scripts reproduce the figures presented in S. Shao, P.E. Jacob, J. Ding, and V. Tarokh (2017) and its supplementary material. Open access to the manuscript and its supplement is also available here.