Title: | Bayesian Mixture Models with JAGS |
---|---|
Description: | Fits finite mixture models of univariate Gaussian distributions using JAGS within a Bayesian framework. |
Authors: | Bettina Gruen [aut, cre] , Martyn Plummer [ctb] |
Maintainer: | Bettina Gruen <[email protected]> |
License: | GPL-2 | GPL-3 |
Version: | 0.7-6 |
Built: | 2024-11-20 06:24:08 UTC |
Source: | CRAN |
Two different plots are currently provided: a plot of different variables against each other and a plot of the same variable against its values in the other classes.
BMMdiag(object, which = 1:2, variables, ask = interactive(), fct1, fct2, xlim, ylim, auto.layout = TRUE, caption = NULL, main = "", ...)
BMMdiag(object, which = 1:2, variables, ask = interactive(), fct1, fct2, xlim, ylim, auto.layout = TRUE, caption = NULL, main = "", ...)
object |
a |
which |
if only one of the plots is required, specify its number. |
variables |
if |
ask |
prompt user before each page of plots |
fct1 |
string: name of transformation function for variable on x-axis. |
fct2 |
string: name of transformation function for variable on y-axis. |
xlim |
if no range for |
ylim |
if no range for |
auto.layout |
logical: if |
caption |
captions to appear above the plots. |
main |
title to each plot (in addition to the above 'caption'). |
... |
further graphical parameters (see 'plot.xy' and 'par') may also be supplied as arguments. |
The plots help determining which variable will induce a unique labelling when taken for ordering of the segments and indicate if the model is overfitted by specifying too many segments.
Bettina Gruen
Creates the BUGS-model specification and return the values for the initialization, the prior specification and the observations.
BMMmodel(y, k, priors, inits = "initsFS", aprioriWeights = 1, no.empty.classes = FALSE, restrict = "none", ...)
BMMmodel(y, k, priors, inits = "initsFS", aprioriWeights = 1, no.empty.classes = FALSE, restrict = "none", ...)
y |
a numeric vector. |
k |
integer indicating the number of segments. |
priors |
specification of priors by a named list or a
|
inits |
specification of initial values by a named list or string indicating the function to be called. |
aprioriWeights |
specification of prior of the a-priori weights. If
|
no.empty.classes |
logical: should it be prevented that empty classes arise during sampling. |
restrict |
one of "none", "mu", "tau". |
... |
further parameters for the function specified in
|
By default the function initsFS
is called for generating
initial values. Any other function specified by inits
is
assumed to have at least x
, k
and restrict
as
input parameters.
The parameter restrict
indicates if a
location-shift model ("tau"), a scale contaminated model
("mu") or a model where both variables vary over components shall
be fitted.
If the logical no.empty.classes
is TRUE
there
are observations added to the model that the classes are not empty. This
signifies that the likelihood when sampling the class affiliations is
changed thus that any data point which is sampled and is the last one
in its class stays there.
If y
is specified an object of class BMMmodel
is
returned with components:
inits |
named list. |
data |
named list. |
bugs |
text for .bug-file with prefix missing. |
If y
is missing an object of class BMMsetup
is
returned containing the parameter specifications. When JAGScall
is called with this object as model argument, BMMmodel
is
called with y
and the other parameters as input arguments
before calling JAGS.
Bettina Gruen
data("fish", package = "bayesmix") model <- BMMmodel(fish, k = 4, priors = list(kind = "independence", parameter = "priorsFish", hierarchical = "tau"), initialValues = list(S0 = 2)) model
data("fish", package = "bayesmix") model <- BMMmodel(fish, k = 4, priors = list(kind = "independence", parameter = "priorsFish", hierarchical = "tau"), initialValues = list(S0 = 2)) model
Given a JAGSrun
object with model of class BMMmodel
the
aposteriori probabilities are determined. If plot=TRUE
, the
resulting object of class BMMposteriori
is plotted.
BMMposteriori(object, class, caption = NULL, plot = TRUE, auto.layout = TRUE, ...)
BMMposteriori(object, class, caption = NULL, plot = TRUE, auto.layout = TRUE, ...)
object |
a |
class |
a vector of integers indicating for which classes the posterior probabilities shall be plotted. The default is all. |
caption |
captions to appear above the plots. |
plot |
logical indicating if a plot shall be made. |
auto.layout |
logical: if |
... |
further graphical parameters may also be supplied as arguments. |
Given a JAGSrun
object with model of class BMMmodel
the a
posteriori probabilities are calculated for the unique data points
with respect to the components specified by class
.
There is a BMMposteriori
object returned which is a list
including the following components
data |
vector of unique data points. |
post |
a matrix including the posteriori probability of the data points for each class. |
Bettina Gruen
This function enables a comfortable creation of
BMMpriors
objects, which can be used for specifying the priors
of a BMMmodel
.
BMMpriors(specification, y, eps = 10^-16)
BMMpriors(specification, y, eps = 10^-16)
specification |
named list including |
y |
a numeric vector. |
eps |
a numeric value indicating the smallest value for flat priors. |
In specification
kind
can be used for specifying
if an "independent" or a "conditionallyconjugate"
prior shall be used. parameter
can be a named list of values for
the prior variables or a function name (e.g., "priorsUncertain",
"priorsFish", "priorsRaftery"). hierarchical
can be NULL
or "tau" if a hierarchical prior shall be taken for
.
mod
is a named list which provides the possibility
to override the values from parameter
.
There is an object BMMpriors
returned with components
name |
vector indicating which kind of prior is specified and if it is an hierarchical prior and if appropriate with respect to which variable. |
var |
list of variables in the prior and their specified values. |
Bettina Gruen
data("fish", package = "bayesmix") priors <- BMMpriors(y = fish)
data("fish", package = "bayesmix") priors <- BMMpriors(y = fish)
A numeric vector containing 15 observations of differences in heights between pairs of self-fertilized and cross-fertilized plants grown in the same condition.
data(darwin)
data(darwin)
A numeric vector of length 15.
Darwin's data set contains two extremely small values. Therefore, this data set can be used for outlier modelling.
Abraham, B. and G. Box (1978) Linear models and spurious observations. Applied Statistics, 27, 131–8.
data("darwin", package = "bayesmix") ## Estimated sample density plot(density(darwin[[1]]), ylim = c(0, 0.02), main = "Outlier modelling") ss <- seq(-100, 100, by = 1) ## Normal density with estimated mean and sd of whole sample lines(ss, dnorm(ss, mean = mean(darwin[[1]]), sd = sd(darwin[[1]])), col = "red") ## Normal density with estimated mean and sd of sample, where the 2 ## extremely small values are removed lines(ss, dnorm(ss, mean = mean(darwin[-c(1:2),1]), sd = sd(darwin[-c(1:2),1])), col = "green")
data("darwin", package = "bayesmix") ## Estimated sample density plot(density(darwin[[1]]), ylim = c(0, 0.02), main = "Outlier modelling") ss <- seq(-100, 100, by = 1) ## Normal density with estimated mean and sd of whole sample lines(ss, dnorm(ss, mean = mean(darwin[[1]]), sd = sd(darwin[[1]])), col = "red") ## Normal density with estimated mean and sd of sample, where the 2 ## extremely small values are removed lines(ss, dnorm(ss, mean = mean(darwin[-c(1:2),1]), sd = sd(darwin[-c(1:2),1])), col = "green")
A numeric vector containing 256 observations of fish lengths.
data(fish)
data(fish)
A numeric vector of length 256.
This data set can be used for modeling unobserved heterogeneity, as it can be assumed that underlying categories present in the data are the age groups to which the fish belong.
D. M. Titterington, A. F. M. Smith and U.E. Makov (1985) Statistical Analysis of Finite Mixture Distributions. Wiley.
data("fish", package = "bayesmix") ss <- seq(-3, 13, by = 0.01) hist(fish[[1]], 20, freq = FALSE, main = "Fish data") lines(ss, dnorm(ss, mean(fish[[1]]), sd(fish[[1]])), col = "red")
data("fish", package = "bayesmix") ss <- seq(-3, 13, by = 0.01) hist(fish[[1]], 20, freq = FALSE, main = "Fish data") lines(ss, dnorm(ss, mean(fish[[1]]), sd(fish[[1]])), col = "red")
Initial values for nodes are created after the suggestion in Sylvia Fruehwirth-Schnatter's book.
initsFS(x, k, restrict, initialValues = list())
initsFS(x, k, restrict, initialValues = list())
x |
a numeric vector. |
k |
number of segments. |
initialValues |
additional initial values specifications. |
restrict |
one of "none", "mu", "tau". |
The initial values for are determined by the quantiles of
the data, those for
give equal weight on each segment and
those for
are equal for all segments and estimated by the
inverse of the IQR of the data divided by 1.34 and squared.
A list with initial values for the parameter indicated by the name of the respective list element is returned.
Bettina Gruen
Using functionality from package rjags the JAGS model is created and posterior draws are sampled.
JAGScall(model, y, prefix, control, ...)
JAGScall(model, y, prefix, control, ...)
model |
|
y |
a numeric vector. |
prefix |
character: prefix for .bug-file. |
control |
named list or |
... |
additional parameters handed over to |
Returns a "JAGSrun"
object.
Bettina Gruen
The information on the number of burn-in draws and monitored draws is specified. Furthermore, it includes the information which variables shall be monitored and possibly a seed and a random number generator (RNG).
JAGScontrol(variables, n.iter = 1000, thin = 1, burn.in = 0, seed, rng = c("base::Wichmann-Hill", "base::Marsaglia-Multicarry", "base::Super-Duper", "base::Mersenne-Twister"))
JAGScontrol(variables, n.iter = 1000, thin = 1, burn.in = 0, seed, rng = c("base::Wichmann-Hill", "base::Marsaglia-Multicarry", "base::Super-Duper", "base::Mersenne-Twister"))
variables |
names of variables which shall be monitored. |
n.iter |
number of monitored draws. |
thin |
thinning interval. |
burn.in |
number of discarded burn-in draws. |
seed |
integer setting the seed for the RNG. |
rng |
specification of RNG. |
An object of class JAGScontrol
is returned which is a list
containing the specified information.
Bettina Gruen
control <- JAGScontrol(variables = "mu") control
control <- JAGScontrol(variables = "mu") control
Calls jags for MCMC sampling.
JAGSrun(y, prefix = yname, model = BMMmodel(k = 2), control = JAGScontrol(variables = c("mu", "tau", "eta")), tmp = TRUE, cleanup = TRUE, ...)
JAGSrun(y, prefix = yname, model = BMMmodel(k = 2), control = JAGScontrol(variables = c("mu", "tau", "eta")), tmp = TRUE, cleanup = TRUE, ...)
y |
a numeric vector. |
prefix |
character: prefix for .bug-file. |
model |
object of class |
control |
specification of control by a |
tmp |
logical: shall the files be written in a temporary directory. |
cleanup |
logical: shall the created files be removed. |
yname |
a character string with the actual |
... |
further parameters handed over to |
This function is a wrapper calling JAGScall
.
Returns a JAGSrun
object with components
call |
the matched call. |
results |
MCMC results obtained by calling functionality from package rjags. |
model |
a |
variables |
vector containing the names of the monitored variables. |
data |
a numeric vector. |
Bettina Gruen
data("fish", package = "bayesmix") prefix <- "fish" variables <- c("mu","tau","eta") k <- 3 modelFish <- BMMmodel(k = k, priors = list(kind = "independence", parameter = "priorsFish", hierarchical = "tau")) controlFish <- JAGScontrol(variables = c(variables, "S"), n.iter = 100) z1 <- JAGSrun(fish, prefix, model = modelFish, initialValues = list(S0 = 2), control = controlFish, cleanup = TRUE, tmp = FALSE) zSort <- Sort(z1, "mu") BMMposteriori(zSort) data("darwin", package = "bayesmix") prefix <- "darwin" k <- 2 modelDarwin <- BMMmodel(k = k, priors = list(kind = "independence", parameter = "priorsUncertain"), aprioriWeights = c(1, 15), no.empty.classes = TRUE, restrict = "tau") z2 <- JAGSrun(darwin, prefix, model = modelDarwin, control = JAGScontrol(variables = variables, n.iter = 3000, burn.in = 1000), cleanup = TRUE, tmp = FALSE) plot(z2, variables = "mu")
data("fish", package = "bayesmix") prefix <- "fish" variables <- c("mu","tau","eta") k <- 3 modelFish <- BMMmodel(k = k, priors = list(kind = "independence", parameter = "priorsFish", hierarchical = "tau")) controlFish <- JAGScontrol(variables = c(variables, "S"), n.iter = 100) z1 <- JAGSrun(fish, prefix, model = modelFish, initialValues = list(S0 = 2), control = controlFish, cleanup = TRUE, tmp = FALSE) zSort <- Sort(z1, "mu") BMMposteriori(zSort) data("darwin", package = "bayesmix") prefix <- "darwin" k <- 2 modelDarwin <- BMMmodel(k = k, priors = list(kind = "independence", parameter = "priorsUncertain"), aprioriWeights = c(1, 15), no.empty.classes = TRUE, restrict = "tau") z2 <- JAGSrun(darwin, prefix, model = modelDarwin, control = JAGScontrol(variables = variables, n.iter = 3000, burn.in = 1000), cleanup = TRUE, tmp = FALSE) plot(z2, variables = "mu")
Plot method for object of class BMMposteriori
, typically
called by BMMposteriori
.
## S3 method for class 'BMMposteriori' plot(x, caption, main = "", ...)
## S3 method for class 'BMMposteriori' plot(x, caption, main = "", ...)
x |
a |
caption |
captions to appear above the plots. |
main |
title to each plot-in addition to the above |
... |
further graphical parameters may also be supplied as arguments. |
This function is called by BMMposteriori
if plot =
TRUE
.
Bettina Gruen
Plots mcmc
chains of a JAGSrun
object.
## S3 method for class 'JAGSrun' plot(x, variables = NULL, trace = TRUE, density = TRUE, smooth = TRUE, bwf, num, xlim, auto.layout = TRUE, ask = interactive(), ...)
## S3 method for class 'JAGSrun' plot(x, variables = NULL, trace = TRUE, density = TRUE, smooth = TRUE, bwf, num, xlim, auto.layout = TRUE, ask = interactive(), ...)
x |
a |
variables |
names of variables which shall be plotted. Default
are all names of results except those with a column dimension larger
than the number of classes |
trace |
plot trace of each variable. |
density |
plot density estimate of each variable. |
smooth |
draw a smooth line through trace plots. |
bwf |
bandwidth function for density plots. |
num |
if not all classes of a variable shall be plotted, a subset can be specified. |
xlim |
if not specified, the range of each variable over all classes is taken as default. |
auto.layout |
automatically generate output format. |
ask |
prompt user before each page of plots. |
... |
further arguments for |
Adapted from plot.mcmc
.
Currently only implemented for JAGSrun
objects with model of class
BMMmodel
. Otherwise the default plot method for the results
of the JAGSrun
object is called (plot.mcmc
).
Bettina Gruen
plot.mcmc
, BMMdiag
,
BMMposteriori
Given the data values for the priors are determined.
priorsFish(y, eps = 10^-16) priorsRaftery(y) priorsUncertain(y, eps = 10^-16)
priorsFish(y, eps = 10^-16) priorsRaftery(y) priorsUncertain(y, eps = 10^-16)
y |
a numeric vector. |
eps |
a numeric value indicating the smallest value for flat priors. |
Values for the prior parameter b0
, B0
, nu0
and S0
are determined.
There is a list returned with named components of the prior parameters.
Bettina Gruen
Random permutation of segment labels for each draw in order to get a better estimate of the unrestricted likelihood.
randomPermutation(x)
randomPermutation(x)
x |
a |
The draws are permuted with respect to the different classes
k
.
The input object with permuted results for each draw is returned.
Any variables where there are neither k
different chains nor only one chain observed are dropped.
Bettina Gruen
Ascending sorting of results of JAGSrun
object with
model of class BMMmodel
with respect to a given variable.
Sort(x, by = NULL)
Sort(x, by = NULL)
x |
a |
by |
variable name according to which the segments shall be ordered. |
If by
is not specified, the first variable in the
corresponding vector of the JAGSrun
object is taken.
The input object with results sorted in ascending order according
to the variable given in by
is returned.
If there arise problems, the original object is returned with a warning.
Bettina Gruen