Performs probability mass function estimation under nonparametric mixture of Poisson kernels.
npmp(y, k, nrep, nb, alpha=1, theta=alpha, sigma=0,
mixing_hyperprior= FALSE, basemeasure_hyperprior = FALSE, mixing_type="DP",
algo="slice", prior="gamma", a, b, a_a=1, b_a=1, lb=NULL, ub=NULL,
print = 1, ndisplay = nrep/4, plot.it = FALSE, pdfwrite = FALSE, ... )dpmpoiss(y, k, nrep, nb, alpha = 1, a, b, lb = NULL, ub = NULL, print = 1,
ndisplay = nrep/4, plot.it = FALSE, pdfwrite = FALSE, ...)
Name of the model
Name of the mixing prior
Quantities about MCMC sampling
MCMC chains of the parameters
A list containing several quantities related to the probability mass function (emprical pmf, posterior mean pmf and pointwise 95% credible intervals) computed for the values from lb to ub
A list containing the posterior mean of the cluster specific parameters (be careful of label-switching problems)
A list containing posterior quantities related to the clustering structure of the data
Vector of count data
Truncation level for the number of cluster in the mixture. Default is length(ydis).
Number of MCMC iterations
Number of burn-in iteration in the MCMC to discard
Value of the precision parameter of the Dirichlet process prior
Value of the strength parameter of the Two-parameters-Poisson-Dirichlet process prior
Value of the discount parameter of the Two-parameters-Poisson-Dirichlet process prior
Logical. If TRUE alpha is random with gamma hyperprior
Logical. If TRUE also the parameters of the base measure are random, see details below.
Type of mixing distribution. Default is "DP" for Dirichlet process but also "2PD" for Two-parameters-Poisson-Dirichlet process is allowed.
Type of algorithm. Current choices are: slice sampler (algo="slice") or polya-urn-type sampler (algo="polya-urn").
String for the base measure prior. Default is "gamma" for lambda ~ Gamma(a,b). The other choice is "normal" for exp(lambda) ~ N(a,b)
Shape (mean) hyperparameter for the gamma (normal) base measure
Scale (sd) hyperparameter for the gamma (normal) prior
Shape hyperparameter for alpha
Scale hyperparameter for alpha
Scalar integer. Lower bound for the argument of the pmf. Default is max(0,min(ydis)-10).
Scalar integer. Upper bound for the argument of the pmf. Default is max(ydis)+10.
Vector of integers (from 1 to 5) indicating whether to print each step of the Gibbs sampler. Specifically, 1 for current iteration, 2 for the DP cluster allocation, 3 for the posterior parameters of the mixture components, 4 for the precision of the DP, 5 for the posterior pmf.
Scalar integer. It gives the number of iterations to be displayed on screen (the function reports on the screen when every ndisplay iterations have been carried out)
Logical, default FALSE. If TRUE a plot with empirical and estimated posterior probability mass functions is plotted.
Logical, default FALSE. If TRUE a pdf file is written in the current working directory. Traceplots and other posterior quantities are drown.
Additional arguments (for future implemetantions).
R code and porting by A. Canale, C code by A. Canale with minor contributions by N. Lunardon.
The function npmp performs probability mass function estimation under nonparametric mixture of Poisson kernels, i.e.
$$y_i \mid \lambda_i \sim \mbox{Poi}(\lambda_i), i=1, \dots, n$$
$$\lambda_i \mid G \sim G$$
$$G \sim \Pi(P_0),$$
where \(\Pi\) is a nonparametric prior (Dirichlet process or Two-parameters-Poisson-Dirichlet process) with base measure \(P_0\). The function dppoiss is a wrapper to npmp with mixing_type="DP" for back portability with version 1.0 of the package. The main part of the code is written in C language to gain computational speed. Plots and posterior summaries are in plain R code. From version 2.0 on, the blocked Gibbs sampler has been removed in place of slice samper (Kalli et al., 2011) and polya-urn sampler. Two different base measures \(P_0\) are implemented: prior="gamma" and prior="normal" for
$$
\lambda_h \sim \mbox{Gamma}(a,b), \log(\lambda_h) \sim N(a,b), h=1,\dots
$$
respectively.
Canale, A. and Dunson, D. B. (2011), "Bayesian Kernel Mixtures for Counts", Journal of American Statistical Association, 106, 1528-1539.
Kalli, M., Griffin, J., and Walker, S. (2011), "Slice sampling mixture models," Statistics and Computing, 21, 93-105.
rmg
# \donttest{
data(ethylene)
y <- tapply(ethylene$impl,FUN=mean,INDEX=ethylene$id)
z <- tapply(ethylene$dose,FUN=mean,INDEX=ethylene$id)
# Estimate the pmf of the number of implants in the control group
y0 <- y[z==0]
pmf.control = dpmpoiss(y0, k=20, nrep=11000, nb=1000, alpha=1, a=1, b=1,
lb=5, ub=24, plot.it=TRUE)
# }
Run the code above in your browser using DataLab