estim.gsm_theta {GSM}

Estimation of a Gamma Shape Mixture Model (GSM)

Description

This function provides the inferential algorithm to estimate a mixture of gamma distributions in which the mixing occurs over the shape parameter. It implements the standard approach for the GSM model, as discussed in Venturini et al. (2008).

Usage

   estim.gsm_theta(y, J, G = 100, M = 600, a, b, alpha, init = list(rep(1 / J, J),
   		J / max(y), rep(1, N)))

Arguments

y

vector of data.

J

number of mixture components.

G

number of points where to evaluate the GSM density.

M

number of MCMC runs.

a

hyperparameter of the rate parameter prior distribution.

b

hyperparameter of the rate parameter prior distribution.

alpha

hyperparameter of the mixture's weights prior distribution.

init

initialization values.

Details

Suggestions on how to choose J, a and b are provided in Venturini et al. (2008). In that work the alpha vector is always set at (1/J,...,1/J), but here one is free to choose the value of the generic element of alpha.

Value

estim.gsm_theta returns an object of class "gsm", which is a list with the following components:

fdens

matrix containing the posterior draws for the mixture's density.

theta

vector containing the posterior draws for the mixture's rate parameter.

weight

matrix containing the posterior draws for the mixture's weights.

label

matrix containing the posterior draws for the mixture's labels.

data

vector of data.

Author(s)

Sergio Venturini [e-mail]

References

Venturini, S., Dominici, F. and Parmigiani, G. (2008), "Gamma shape mixtures for heavy-tailed distributions". Annals of Applied Statistics, Volume 2, Number 2, 756–776. http://projecteuclid.org/euclid.aoas/1215118537

See Also

estim.gsm, summary-methods, plot-methods.

Examples

## Not run: 
set.seed(2040)
y <- rgsm(500, c(.1, .3, .4, .2), 1)
burnin <- 100
mcmcsim <- 500
J <- 250
gsm.out <- estim.gsm_theta(y, J, 300, burnin + mcmcsim, 6500, 340, 1/J)
summary(gsm.out, plot = TRUE, start = (burnin + 1))
plot(gsm.out, ndens = 0, nbin = 20, histogram = TRUE, start = (burnin + 1))
## End(Not run)