Title: | Ferguson-Klass Type Algorithm for Posterior Normalized Random Measures |
---|---|
Description: | Bayesian nonparametric density estimation modeling mixtures by a Ferguson-Klass type algorithm for posterior normalized random measures. |
Authors: | Julyan Arbel [ctb], Ernesto Barrios [aut], Guillaume Kon Kam King [aut, cre], Antonio Lijoi [aut], Luis E. Nieto-Barajas [aut], Igor Prünster [aut] |
Maintainer: | Guillaume Kon Kam King <[email protected]> |
License: | GPL (>= 2) |
Version: | 2023.3.8 |
Built: | 2024-12-18 06:33:23 UTC |
Source: | CRAN |
This package performs Bayesian nonparametric density estimation for exact and censored data via a normalized random measure mixture model. The package allows the user to specify the mixture kernel, the mixing normalized measure and the choice of performing fully nonparametric mixtures on locations and scales, or semiparametric mixtures on locations only with common scale parameter. Options for the kernels are: two kernels with support in the real line (gaussian and double exponential), two more kernels in the positive line (gamma and lognormal) and one with bounded support (beta). The options for the normalized random measures are members of the class of normalized generalized gamma, which include the Dirichlet process, the normalized inverse gaussian process and the normalized stable process. The type of censored data handled by the package is right, left and interval.
Package: | BNPdensity |
Type: | Package |
Version: | 2016.10 |
Date: | 2016-10-14 |
License: | GPL version 2 or later |
LazyLoad: | yes |
The package includes four main functions: MixNRMI1, MixNRMI2, MixNRMI1cens and MixNRMI2cens which implement semiparametric and fully nonparametric mixtures for exact data, and semiparametric and fully nonparametric mixtures for censored data respectively. Additionally, the package includes several other functions required for sampling from conditional distributions in the MCMC implementation. These functions are intended for internal use only.
Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I.; Contributor: Guillaume Kon Kam King.; Maintainer: Ernesto Barrios <ebarrios at itam.mx>
Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I. (2013). Modeling with Normalized Random Measure Mixture Models. Statistical Science. Vol. 28, No. 3, 313-334.
Kon Kam King, G., Arbel, J. and Prünster, I. (2016). Species Sensitivity Distribution revisited: a Bayesian nonparametric approach. In preparation.
MixNRMI1
, MixNRMI2
,
MixNRMI1cens
, MixNRMI2cens
example(MixNRMI1) example(MixNRMI2) example(MixNRMI1cens) example(MixNRMI2cens)
example(MixNRMI1) example(MixNRMI2) example(MixNRMI1cens) example(MixNRMI2cens)
Concerns an acidity index measured in a sample of 155 lakes in north-central Wisconsin.
A real vector with 155 observations.
Crawford, S. L., DeGroot, M. H., Kadane, J. B. and Small, M. J. (1992). Modeling lake chemistry distributions: approximate Bayesian methods for estimating a finite mixture model. Technometrics, 34, 441-453.
data(acidity) hist(acidity)
data(acidity) hist(acidity)
This is a helper function for use in Reduce() over a list of vectors
add(x, y)
add(x, y)
x |
first argument of the sum |
y |
second argument of the sum |
x + y
Convert the output of multMixNRMI into a coda mcmc object
## S3 method for class 'multNRMI' as.mcmc(x, ..., thinning_to = 1000, ncores = parallel::detectCores())
## S3 method for class 'multNRMI' as.mcmc(x, ..., thinning_to = 1000, ncores = parallel::detectCores())
x |
Output of multMixNRMI. |
... |
Further arguments to be passed to specific methods |
thinning_to |
Final length of the chain after thinning. |
ncores |
Specify the number of cores to use in the conversion |
a coda::mcmc object
data(acidity) out <- multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2) coda::as.mcmc(out, ncores = 2)
data(acidity) out <- multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2) coda::as.mcmc(out, ncores = 2)
The function Rmpfr::asNumeric prints the following warning: In asMethod(object) : coercing "mpfr1" via "mpfr" (inefficient). It is not clear how to avoid it nor how to silence it, hence this function. A cleaner solution may be available at: https://stackoverflow.com/questions/4948361/how-do-i-save-warnings-and-errors-as-output-from-a-function/4952908#4952908
asNumeric_no_warning(x)
asNumeric_no_warning(x)
x |
An object of class Rmpfr::mpfr1 |
a "numeric" number
Comment on the NRMI process depending on the value of the parameters
comment_on_NRMI_type(NRMI_param = list(Alpha = 1, Kappa = 0, Gamma = 0.4))
comment_on_NRMI_type(NRMI_param = list(Alpha = 1, Kappa = 0, Gamma = 0.4))
NRMI_param |
A named list of the form list("Alpha" = 1, "Kappa" = 0, "Gamma" = 0.4) |
A string containing a comment on the NRMI process
BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0, "Gamma" = 0.4)) BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0.1, "Gamma" = 0.4)) BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0.1, "Gamma" = 0.5))
BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0, "Gamma" = 0.4)) BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0.1, "Gamma" = 0.4)) BNPdensity:::comment_on_NRMI_type(list("Alpha" = 1, "Kappa" = 0.1, "Gamma" = 0.5))
Summarizes the posterior on all possible clusterings by an optimal clustering where optimality is defined as minimizing the posterior expectation of a specific loss function, the Variation of Information or Binder's loss function. Computation can be lengthy for large datasets, because of the large size of the space of all clusterings.
compute_optimal_clustering(fit, loss_type = "VI")
compute_optimal_clustering(fit, loss_type = "VI")
fit |
The fitted object, obtained from one of the MixNRMIx functions |
loss_type |
Defines the loss function to be used in the expected posterior loss minimization. Can be one of "VI" (Variation of Information), "B" (Binder's loss), "NVI" (Normalized Variation of Information) or "NID" (Normalized Information Distance). Defaults to "VI". |
A vector of integers with the same size as the data, indicating the allocation of each data point.
This function creates an real grid then rounds it. If the grid is fine enough, there is a risk that rounding ties, i.e. iteration which are kept twice. To avoid this, if the total number of iterations is smaller than twice the number of iterations desired after thinning, the chain is not thinned.
compute_thinning_grid(Nit, thinning_to = 10)
compute_thinning_grid(Nit, thinning_to = 10)
Nit |
Length of the MCMC chain |
thinning_to |
Desired number of iterations after thinning. |
an integer vector of the MCMC iterations retained.
Convert the output of multMixNRMI into a coda mcmc object
convert_to_mcmc(fitlist, thinning_to = 1000, ncores = parallel::detectCores())
convert_to_mcmc(fitlist, thinning_to = 1000, ncores = parallel::detectCores())
fitlist |
Output of multMixNRMI. |
thinning_to |
Final length of the chain after thinning. |
ncores |
Specify the number of cores to use in the conversion |
a coda::mcmc object
This function assumes that all chains have the same size. To allow for different chain sizes, care should be paid to proper weighting.
## S3 method for class 'multNRMI' cpo(object, ...)
## S3 method for class 'multNRMI' cpo(object, ...)
object |
A fit obtained through from the functions MixNRMI1/MixNRMI1cens |
... |
Further arguments to be passed to generic function, ignored at the moment |
A vector of Conditional Predictive Ordinates (CPOs)
data(acidity) out <- multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2) cpo(out)
data(acidity) out <- multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2) cpo(out)
Extract the Conditional Predictive Ordinates (CPOs) from a fitted object
## S3 method for class 'NRMI1' cpo(object, ...)
## S3 method for class 'NRMI1' cpo(object, ...)
object |
A fit obtained through from the functions MixNRMI1/MixNRMI1cens |
... |
Further arguments to be passed to generic function, ignored at the moment |
A vector of Conditional Predictive Ordinates (CPOs)
data(acidity) out <- MixNRMI1(acidity, Nit = 50) cpo(out)
data(acidity) out <- MixNRMI1(acidity, Nit = 50) cpo(out)
Extract the Conditional Predictive Ordinates (CPOs) from a fitted object
## S3 method for class 'NRMI2' cpo(object, ...)
## S3 method for class 'NRMI2' cpo(object, ...)
object |
A fit obtained through from the function MixNRMI2/MixNRMI2cens |
... |
Further arguments to be passed to generic function, ignored at the moment |
A vector of Conditional Predictive Ordinates (CPOs)
data(acidity) out <- MixNRMI2(acidity, Nit = 50) cpo(out)
data(acidity) out <- MixNRMI2(acidity, Nit = 50) cpo(out)
Convert distribution names to indices
dist_name_k_index_converter(distname)
dist_name_k_index_converter(distname)
distname |
a character representing the distribution name. Allowed names are "normal", "gamma", "beta", "exponential", "double exponential", "lognormal", "half-Cauchy", "half-normal", "half-student", "uniform" and "truncated normal", or their common abbreviations "norm", "exp", "lnorm", "halfcauchy", "halfnorm", "halft" and "unif". |
an index describing the distribution. 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal, 6 = Half-Cauchy, 7 = Half-normal, 8 = Half-Student, 9 = Uniform, 10 = Truncated normal
Computes the density.
dt_(x, df, mean, sd)
dt_(x, df, mean, sd)
x |
Numeric vector. Data set to which the density is evaluated. |
df |
Numeric constant. Degrees of freedom (> 0, maybe non-integer) |
mean |
Numeric constant. Location parameter. |
sd |
Positive numeric constant. Scale parameter. ## The function is currently defined as function(x, df, mean, sd) dt((x - mean) / sd, df, ncp = 0) / sd |
For internal use
Concerns the distribution of enzymatic activity in the blood, for an enzyme involved in the metabolism of carcinogenetic substances, among a group of 245 unrelated individuals.
A data frame with 244 observations on the following variable:
A numeric vector.
Bechtel, Y. C., Bonaiti-Pellie, C., Poisson, N., Magnette, J. and Bechtel, P.R. (1993). A population and family study of N-acetyltransferase using caffeine urinary metabolites. Clin. Pharm. Therp., 54, 134-141.
data(enzyme) hist(enzyme)
data(enzyme) hist(enzyme)
This object contains the output when setting set.seed(150520) and running the function Enzyme1.out <- MixNRMI1(enzyme, Alpha = 1, Kappa = 0.007, Gama = 0.5, distr.k = "gamma", distr.p0 = "gamma", asigma = 1, bsigma = 1, Meps = 0.005, Nit = 5000, Pbi = 0.2)
See function MixNRMI1
data(Enzyme1.out)
data(Enzyme1.out)
This object contains the output when setting set.seed(150520) and running the function Enzyme2.out <- MixNRMI2(enzyme, Alpha = 1, Kappa = 0.007, Gama = 0.5, distr.k = "gamma", distr.py0 = "gamma", distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps = 0.005, Nit = 5000, Pbi = 0.2) See function MixNRMI2
data(Enzyme2.out)
data(Enzyme2.out)
Computes the expected number of components for a Dirichlet process.
expected_number_of_components_Dirichlet( n, Alpha, ntrunc = NULL, silence = TRUE )
expected_number_of_components_Dirichlet( n, Alpha, ntrunc = NULL, silence = TRUE )
n |
Number of data points |
Alpha |
Numeric constant. Total mass of the centering measure. |
ntrunc |
Level of truncation when computing the expectation. Defaults to n. If greater than n, it is fixed to n. |
silence |
Boolean. Whether to print the current calculation step for the Stable process, as the function can be long |
A real value which approximates the expected number of components
Reference: P. De Blasi, S. Favaro, A. Lijoi, R. H. Mena, I. Prünster, and M. Ruggiero, “Are Gibbs-type priors the most natural generalization of the Dirichlet process?,” IEEE Trans. Pattern Anal. Mach. Intell., vol. 37, no. 2, pp. 212–229, 2015.
expected_number_of_components_Dirichlet(100, 1.2)
expected_number_of_components_Dirichlet(100, 1.2)
Computes the expected number of components for a stable process.
expected_number_of_components_stable(n, Gama, ntrunc = NULL)
expected_number_of_components_stable(n, Gama, ntrunc = NULL)
n |
Number of data points |
Gama |
Numeric constant. 0 <= Gama <=1. |
ntrunc |
Level of truncation when computing the expectation. Defaults to n. If greater than n, it is fixed to n. |
A real value of type mpfr1 which approximates the expected number of components
In spite of the high precision arithmetic packages used for in function, it can be numerically unstable for small values of Gama. This is because evaluating a sum with alternated signs, in the generalized factorial coefficients, is tricky. Reference: P. De Blasi, S. Favaro, A. Lijoi, R. H. Mena, I. Prünster, and M. Ruggiero, “Are gibbs-type priors the most natural generalization of the Dirichlet process?,” IEEE Trans. Pattern Anal. Mach. Intell., vol. 37, no. 2, pp. 212–229, 2015.
expected_number_of_components_stable(100, 0.8)
expected_number_of_components_stable(100, 0.8)
Repeat the common scale parameter of a semiparametric model to match the dimension of the location parameters.
fill_sigmas(semiparametric_fit)
fill_sigmas(semiparametric_fit)
semiparametric_fit |
The result of the fit, obtained through the function MixNRMI1. |
an adequate list of vectors of sigmas
Velocities of 82 galaxies diverging from our own galaxy.
A data frame with 82 observations on the following variable:
A numeric vector.
Roeder, K. (1990) "Density estimation with confidence sets exemplified by superclusters and voids in the galaxies". Journal of the American Statistical Association. 85, 617-624.
data(galaxy) hist(galaxy)
data(galaxy) hist(galaxy)
This object contains the output when setting set.seed(150520) and running the function MixNRMI1(galaxy, Alpha = 1, Kappa = 0.015, Gama = 0.5, distr.k = "normal", distr.p0 = "gamma", asigma = 1, bsigma = 1, delta = 7, Meps = 0.005, Nit = 5000, Pbi = 0.2)
See function MixNRMI1.
data(Galaxy1.out)
data(Galaxy1.out)
This object contains the output when setting set.seed(150520) and running the function Enzyme2.out <- MixNRMI2(x, Alpha = 1, Kappa = 0.007, Gama = 0.5, distr.k = "gamma", distr.py0 = "gamma", distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps = 0.005, Nit = 5000, Pbi = 0.2)
See function MixNRMI2.
data(Galaxy2.out)
data(Galaxy2.out)
This function is used in the print methods for MixNRMI1, MixNRMI2, MixNRMI1cens, MixNRMI2cens, and all the multMixNRMIx versions
give_kernel_name(distr.k)
give_kernel_name(distr.k)
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
A character with the name of the distribution used as the kernel
BNPdensity:::give_kernel_name(4)
BNPdensity:::give_kernel_name(4)
Plot Goodness of fits graphical checks for censored data
GOFplots(fit, qq_plot = FALSE, thinning_to = 500)
GOFplots(fit, qq_plot = FALSE, thinning_to = 500)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2, MixMRMI1cens or MixMRMI2cens |
qq_plot |
Whether to compute the QQ-plot |
thinning_to |
How many iterations to compute the mean posterior quantiles |
A density plot, a cumulative density plot with the Turnbull cumulative distribution, a percentile-percentile plot, and potentially a quantile-quantile plot.
set.seed(150520) data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) GOFplots(out)
set.seed(150520) data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) GOFplots(out)
Plot Goodness of fits graphical checks for censored data
GOFplots_censored(fit, qq_plot = FALSE, thinning_to = 500)
GOFplots_censored(fit, qq_plot = FALSE, thinning_to = 500)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2, MixMRMI1cens or MixMRMI2cens |
qq_plot |
Whether to compute the QQ-plot |
thinning_to |
How many iterations to compute the mean posterior quantiles |
A density plot, a cumulative density plot with the Turnbull cumulative distribution, and a percentile-percentile plot.
set.seed(150520) data(salinty) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) BNPdensity:::GOFplots_censored(out)
set.seed(150520) data(salinty) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) BNPdensity:::GOFplots_censored(out)
Plot Goodness of fits graphical checks for non censored data
GOFplots_noncensored(fit, qq_plot = FALSE, thinning_to = 500)
GOFplots_noncensored(fit, qq_plot = FALSE, thinning_to = 500)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2, MixMRMI1cens or MixMRMI2cens |
qq_plot |
Whether to compute the QQ-plot |
thinning_to |
How many iterations to compute the mean posterior quantiles |
A density plot with histogram, a cumulative density plot with the empirical cumulative distribution, and a percentile-percentile plot.
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::GOFplots_noncensored(out)
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::GOFplots_noncensored(out)
Create a plotting grid from censored or non-censored data.
grid_from_data(data, npoints = 100)
grid_from_data(data, npoints = 100)
data |
Input data from which to compute the grid. |
npoints |
Number of points on the grid. |
a vector containing the plotting grid
Create a plotting grid from censored data.
grid_from_data_censored(data, npoints = 100)
grid_from_data_censored(data, npoints = 100)
data |
Censored input data from which to compute the grid. |
npoints |
Number of points on the grid. |
a vector containing the plotting grid
Create a plotting grid from non-censored data.
grid_from_data_noncensored(data, npoints = 100)
grid_from_data_noncensored(data, npoints = 100)
data |
Non-censored input data from which to compute the grid. |
npoints |
Number of points on the grid. |
a vector containing the plotting grid
Test if the data is censored
is_censored(dat)
is_censored(dat)
dat |
The dataset to be tested |
TRUE if the data is censored
data(salinity) BNPdensity:::is_censored(salinity)
data(salinity) BNPdensity:::is_censored(salinity)
Tests if a fit is a semi parametric or nonparametric model.
is_semiparametric(fit)
is_semiparametric(fit)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2. |
TRUE if the fit is a semiparametric model
set.seed(150520) data(acidity) x <- enzyme out <- MixNRMI1(enzyme, extras = TRUE, Nit = 10) BNPdensity:::is_semiparametric(out)
set.seed(150520) data(acidity) x <- enzyme out <- MixNRMI1(enzyme, extras = TRUE, Nit = 10) BNPdensity:::is_semiparametric(out)
Bayesian nonparametric estimation based on normalized measures driven mixtures for locations.
MixNRMI1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = 1, asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
MixNRMI1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = 1, asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.p0 |
The distribution name for the centering measure. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure: 1 = Normal; 2 = Gamma; 3 = Beta. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior
on the standard deviation of the mixture kernel |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior
on the standard deviation of the mixture kernel |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling sigma. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
This generic function fits a normalized random measure (NRMI) mixture model for density estimation (James et al. 2009). Specifically, the model assumes a normalized generalized gamma (NGG) prior for the locations (means) of the mixture kernel and a parametric prior for the common smoothing parameter sigma, leading to a semiparametric mixture model.
The details of the model are:
where
's are the observed data,
's are latent (location)
variables,
sigma
is the smoothing parameter, k
is a parametric
kernel parameterized in terms of mean and standard deviation, (Alpha,
Kappa, Gama; P_0)
are the parameters of the NGG prior with P_0
being
the centering measure whose parameters are assigned vague hyper prior
distributions, and (asigma,bsigma)
are the hyper-parameters of the
gamma prior on the smoothing parameter sigma
. In particular:
NGG(Alpha, 1, 0; P_0)
defines a Dirichlet process; NGG(1,
Kappa, 1/2; P_0)
defines a Normalized inverse Gaussian process; and
NGG(1, 0, Gama; P_0)
defines a normalized stable process.
The evaluation grid ranges from min(x) - epsilon
to max(x) +
epsilon
. By default epsilon=sd(x)/4
.
The function returns a MixNRMI1 object. It is based on a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
cpo |
Numeric vector of |
R |
Numeric vector of
|
S |
Numeric vector of |
U |
Numeric vector of
|
Allocs |
List of |
means |
List of |
weights |
List of
|
Js |
List of |
Nm |
Integer constant. Number of jumps of the continuous component of the unnormalized process. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
procTime |
Numeric vector with execution time provided by
|
distr.k |
Integer corresponding to the kernel chosen for the mixture |
data |
Data used for the fit |
NRMI_params |
A named list with the parameters of the NRMI process |
The function is computing intensive. Be patient.
Barrios, E., Kon Kam King, G., Lijoi, A., Nieto-Barajas, L.E. and Prüenster, I.
1.- Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I. (2013). Modeling with Normalized Random Measure Mixture Models. Statistical Science. Vol. 28, No. 3, 313-334.
2.- James, L.F., Lijoi, A. and Prünster, I. (2009). Posterior analysis for normalized random measure with independent increments. Scand. J. Statist 36, 76-97.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
### Example 1 ## Not run: # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI1(x) # Plotting density estimate + 95% credible interval plot(out) ### Example 2 set.seed(150520) data(enzyme) x <- enzyme Enzyme1.out <- MixNRMI1(x, Alpha = 1, Kappa = 0.007, Gama = 0.5, distr.k = "gamma", distr.p0 = "gamma", asigma = 1, bsigma = 1, Meps=0.005, Nit = 5000, Pbi = 0.2) attach(Enzyme1.out) # Plotting density estimate + 95% credible interval plot(Enzyme1.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting sigma par(mfrow = c(2, 1)) plot(S, type = "l", main = "Trace of sigma") hist(S, nclass = 20, probability = TRUE, main = "Histogram of sigma") # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach() ## End(Not run) ### Example 3 ## Do not run # set.seed(150520) # data(galaxy) # x <- galaxy # Galaxy1.out <- MixNRMI1(x, Alpha = 1, Kappa = 0.015, Gama = 0.5, # distr.k = "normal", distr.p0 = "gamma", # asigma = 1, bsigma = 1, delta = 7, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(galaxy) x <- galaxy data(Galaxy1.out) attach(Galaxy1.out) # Plotting density estimate + 95% credible interval plot(Galaxy1.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting sigma par(mfrow = c(2, 1)) plot(S, type = "l", main = "Trace of sigma") hist(S, nclass = 20, probability = TRUE, main = "Histogram of sigma") # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach()
### Example 1 ## Not run: # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI1(x) # Plotting density estimate + 95% credible interval plot(out) ### Example 2 set.seed(150520) data(enzyme) x <- enzyme Enzyme1.out <- MixNRMI1(x, Alpha = 1, Kappa = 0.007, Gama = 0.5, distr.k = "gamma", distr.p0 = "gamma", asigma = 1, bsigma = 1, Meps=0.005, Nit = 5000, Pbi = 0.2) attach(Enzyme1.out) # Plotting density estimate + 95% credible interval plot(Enzyme1.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting sigma par(mfrow = c(2, 1)) plot(S, type = "l", main = "Trace of sigma") hist(S, nclass = 20, probability = TRUE, main = "Histogram of sigma") # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach() ## End(Not run) ### Example 3 ## Do not run # set.seed(150520) # data(galaxy) # x <- galaxy # Galaxy1.out <- MixNRMI1(x, Alpha = 1, Kappa = 0.015, Gama = 0.5, # distr.k = "normal", distr.p0 = "gamma", # asigma = 1, bsigma = 1, delta = 7, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(galaxy) x <- galaxy data(Galaxy1.out) attach(Galaxy1.out) # Plotting density estimate + 95% credible interval plot(Galaxy1.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting sigma par(mfrow = c(2, 1)) plot(S, type = "l", main = "Trace of sigma") hist(S, nclass = 20, probability = TRUE, main = "Histogram of sigma") # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach()
Bayesian nonparametric estimation based on normalized measures driven mixtures for locations.
MixNRMI1cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
MixNRMI1cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
xleft |
Numeric vector. Lower limit of interval censoring. For exact data the same as xright |
xright |
Numeric vector. Upper limit of interval censoring. For exact data the same as xleft. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.p0 |
The distribution name for the centering measure. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure: 1 = Normal; 2 = Gamma; 3 = Beta. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior
on the standard deviation of the mixture kernel |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior
on the standard deviation of the mixture kernel |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling sigma. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
This generic function fits a normalized random measure (NRMI) mixture model for density estimation (James et al. 2009) with censored data. Specifically, the model assumes a normalized generalized gamma (NGG) prior for the locations (means) of the mixture kernel and a parametric prior for the common smoothing parameter sigma, leading to a semiparametric mixture model.
This function coincides with MixNRMI1
when the lower (xleft)
and upper (xright) censoring limits correspond to the same exact value.
The details of the model are:
where
's are the observed data,
's are latent (location)
variables,
sigma
is the smoothing parameter, k
is a parametric
kernel parameterized in terms of mean and standard deviation, (Alpha,
Kappa, Gama; P_0)
are the parameters of the NGG prior with P_0
being
the centering measure whose parameters are assigned vague hyper prior
distributions, and (asigma,bsigma)
are the hyper-parameters of the
gamma prior on the smoothing parameter sigma
. In particular:
NGG(Alpha, 1, 0; P_0)
defines a Dirichlet process; NGG(1,
Kappa, 1/2; P_0)
defines a Normalized inverse Gaussian process; and
NGG(1, 0, Gama; P_0)
defines a normalized stable process.
The evaluation grid ranges from min(x) - epsilon
to max(x) +
epsilon
. By default epsilon=sd(x)/4
.
The function returns a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
cpo |
Numeric vector of |
R |
Numeric vector of
|
S |
Numeric vector of |
U |
Numeric vector of
|
Allocs |
List of |
means |
List of |
weights |
List of
|
Js |
List of |
Nm |
Integer constant. Number of jumps of the continuous component of the unnormalized process. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
procTime |
Numeric vector with execution time provided by
|
distr.k |
Integer corresponding to the kernel chosen for the mixture |
data |
Data used for the fit |
NRMI_params |
A named list with the parameters of the NRMI process |
The function is computing intensive. Be patient.
Barrios, E., Kon Kam King, G. and Nieto-Barajas, L.E.
1.- Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I. (2013). Modeling with Normalized Random Measure Mixture Models. Statistical Science. Vol. 28, No. 3, 313-334.
2.- James, L.F., Lijoi, A. and Prünster, I. (2009). Posterior analysis for normalized random measure with independent increments. Scand. J. Statist 36, 76-97.
3.- Kon Kam King, G., Arbel, J. and Prünster, I. (2016). Species Sensitivity Distribution revisited: a Bayesian nonparametric approach. In preparation.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
### Example 1 ## Not run: # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI1cens(x, x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ## Not run: ### Example 2 # Data data(salinity) # Fitting the model under default specifications out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, Nit = 5000) # Plotting density estimate + 95% credible interval attach(out) plot(out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) detach() ## End(Not run)
### Example 1 ## Not run: # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI1cens(x, x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ## Not run: ### Example 2 # Data data(salinity) # Fitting the model under default specifications out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, Nit = 5000) # Plotting density estimate + 95% credible interval attach(out) plot(out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) detach() ## End(Not run)
Bayesian nonparametric estimation based on normalized measures driven mixtures for locations and scales.
MixNRMI2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
MixNRMI2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.py0 |
The distribution name for the centering measure for locations. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure for locations: 1 = Normal; 2 = Gamma; 3 = Beta. |
distr.pz0 |
The distribution name for the centering measure for scales. Allowed names are "gamma", or an integer number identifying the centering measure for
scales: 2 = Gamma. For more options use |
mu.pz0 |
Numeric constant. Prior mean of the centering measure for scales. |
sigma.pz0 |
Numeric constant. Prior standard deviation of the centering measure for scales. |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the scales. |
kappa |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the location parameters. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. If 'adaptive=TRUE', 'delta_U'is the starting value for the adaptation. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, sigmas, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
This generic function fits a normalized random measure (NRMI) mixture model for density estimation (James et al. 2009). Specifically, the model assumes a normalized generalized gamma (NGG) prior for both, locations (means) and standard deviations, of the mixture kernel, leading to a fully nonparametric mixture model.
The details of the model are:
where, 's are the observed data,
's are
bivariate latent (location and scale) vectors,
k
is a parametric
kernel parameterized in terms of mean and standard deviation, (Alpha,
Kappa, Gama; P_0)
are the parameters of the NGG prior with a bivariate
P_0
being the centering measure with independent components, that is,
. The parameters of
P_0(Y)
are assigned
vague hyper prior distributions and (mu.pz0,sigma.pz0)
are the
hyper-parameters of P_0(Z)
. In particular, NGG(Alpha, 1, 0;
P_0)
defines a Dirichlet process; NGG(1, Kappa, 1/2;P_0)
defines a
Normalized inverse Gaussian process; and NGG(1, 0, Gama; P_0)
defines
a normalized stable process. The evaluation grid ranges from min(x) -
epsilon
to max(x) + epsilon
. By default epsilon=sd(x)/4
.
The function returns a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
cpo |
Numeric vector of |
R |
Numeric vector of
|
U |
Numeric vector of |
Allocs |
List of
|
means |
List of |
sigmas |
Numeric vector of
|
weights |
List of |
Js |
List of
|
Nm |
Integer constant. Number of jumps of the continuous component of the unnormalized process. |
delta_Us |
List of
|
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
procTime |
Numeric vector with execution time provided by |
distr.k |
Integer corresponding to the kernel chosen for the mixture |
data |
Data used for the fit |
NRMI_params |
A named list with the parameters of the NRMI process |
The function is computing intensive. Be patient.
Barrios, Kon Kam King, G., E., Lijoi, A., Nieto-Barajas, L.E. and Prüenster, I.
1.- Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I. (2013). Modeling with Normalized Random Measure Mixture Models. Statistical Science. Vol. 28, No. 3, 313-334.
2.- James, L.F., Lijoi, A. and Prünster, I. (2009). Posterior analysis for normalized random measure with independent increments. Scand. J. Statist 36, 76-97.
3.- Arbel, J., Kon Kam King, G., Lijoi, A., Nieto-Barajas, L.E. and Prüenster, I. (2021). BNPdensity: a package for Bayesian Nonparametric density estimation using Normalised Random Measures with Independent Increments.. Australian and New Zealand Journal of Statistics, to appear
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
## Not run: ### Example 1 # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI2(x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ### Example 2 ## Do not run # set.seed(150520) # data(enzyme) # x <- enzyme # Enzyme2.out <- MixNRMI2(x, Alpha = 1, Kappa = 0.007, Gama = 0.5, # distr.k = "gamma", distr.py0 = "gamma", # distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(enzyme) x <- enzyme data(Enzyme2.out) attach(Enzyme2.out) # Plotting density estimate + 95% credible interval plot(Enzyme2.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach() ### Example 3 ## Do not run # set.seed(150520) # data(galaxy) # x <- galaxy # Galaxy2.out <- MixNRMI2(x, Alpha = 1, Kappa = 0.015, Gama = 0.5, # distr.k = "normal", distr.py0 = "gamma", # distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(galaxy) x <- galaxy data(Galaxy2.out) attach(Galaxy2.out) # Plotting density estimate + 95% credible interval plot(Galaxy2.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach()
## Not run: ### Example 1 # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI2(x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ### Example 2 ## Do not run # set.seed(150520) # data(enzyme) # x <- enzyme # Enzyme2.out <- MixNRMI2(x, Alpha = 1, Kappa = 0.007, Gama = 0.5, # distr.k = "gamma", distr.py0 = "gamma", # distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(enzyme) x <- enzyme data(Enzyme2.out) attach(Enzyme2.out) # Plotting density estimate + 95% credible interval plot(Enzyme2.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach() ### Example 3 ## Do not run # set.seed(150520) # data(galaxy) # x <- galaxy # Galaxy2.out <- MixNRMI2(x, Alpha = 1, Kappa = 0.015, Gama = 0.5, # distr.k = "normal", distr.py0 = "gamma", # distr.pz0 = "gamma", mu.pz0 = 1, sigma.pz0 = 1, Meps=0.005, # Nit = 5000, Pbi = 0.2) # The output of this run is already loaded in the package # To show results run the following # Data data(galaxy) x <- galaxy data(Galaxy2.out) attach(Galaxy2.out) # Plotting density estimate + 95% credible interval plot(Galaxy2.out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) # Plotting u par(mfrow = c(2, 1)) plot(U, type = "l", main = "Trace of U") hist(U, nclass = 20, probability = TRUE, main = "Histogram of U") # Plotting cpo par(mfrow = c(2, 1)) plot(cpo, main = "Scatter plot of CPO's") boxplot(cpo, horizontal = TRUE, main = "Boxplot of CPO's") print(paste("Average log(CPO)=", round(mean(log(cpo)), 4))) print(paste("Median log(CPO)=", round(median(log(cpo)), 4))) detach()
Bayesian nonparametric estimation based on normalized measures driven mixtures for locations and scales.
MixNRMI2cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
MixNRMI2cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE )
xleft |
Numeric vector. Lower limit of interval censoring. For exact data the same as xright |
xright |
Numeric vector. Upper limit of interval censoring. For exact data the same as xleft. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.py0 |
The distribution name for the centering measure for locations. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure for locations: 1 = Normal; 2 = Gamma; 3 = Beta. |
distr.pz0 |
The distribution name for the centering measure for scales. Allowed names are "gamma", "lognormal", "half-Cauchy", "half-normal", "half-student", "uniform" and "truncated normal", or their common abbreviations "norm", "exp", "lnorm", "halfcauchy", "halfnorm", "halft" and "unif", or an integer number identifying the centering measure for scales: 2 = Gamma, 5 = Lognormal, 6 = Half Cauchy, 7 = Half Normal, 8 = Half Student-t, 9 = Uniform, 10 = Truncated Normal. |
mu.pz0 |
Numeric constant. Prior mean of the centering measure for scales. |
sigma.pz0 |
Numeric constant. Prior standard deviation of the centering measure for scales. |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the scales. |
kappa |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the location parameters. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. If 'adaptive=TRUE', 'delta_U'is the starting value for the adaptation. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, sigmas, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
This generic function fits a normalized random measure (NRMI) mixture model for density estimation (James et al. 2009). Specifically, the model assumes a normalized generalized gamma (NGG) prior for both, locations (means) and standard deviations, of the mixture kernel, leading to a fully nonparametric mixture model.
The details of the model are:
where, 's are the observed data,
's are
bivariate latent (location and scale) vectors,
k
is a parametric
kernel parameterized in terms of mean and standard deviation, (Alpha,
Kappa, Gama; P_0)
are the parameters of the NGG prior with a bivariate
P_0
being the centering measure with independent components, that is,
. The parameters of
P_0(Y)
are assigned
vague hyper prior distributions and (mu.pz0,sigma.pz0)
are the
hyper-parameters of P_0(Z)
. In particular, NGG(Alpha, 1, 0;
P_0)
defines a Dirichlet process; NGG(1, Kappa, 1/2;P_0)
defines a
Normalized inverse Gaussian process; and NGG(1, 0, Gama; P_0)
defines
a normalized stable process. The evaluation grid ranges from min(x) -
epsilon
to max(x) + epsilon
. By default epsilon=sd(x)/4
.
The function returns a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
cpo |
Numeric vector of |
R |
Numeric vector of
|
U |
Numeric vector of |
Allocs |
List of
|
means |
List of |
sigmas |
Numeric vector of
|
weights |
List of |
Js |
List of
|
Nm |
Integer constant. Number of jumps of the continuous component of the unnormalized process. |
delta_Us |
List of
|
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric
constant. Burn-in period proportion of |
procTime |
Numeric
vector with execution time provided by |
distr.k |
Integer corresponding to the kernel chosen for the mixture |
data |
Data used for the fit |
NRMI_params |
A named list with the parameters of the NRMI process |
The function is computing intensive. Be patient.
Barrios, E., Kon Kam King, G. and Nieto-Barajas, L.E.
1.- Barrios, E., Lijoi, A., Nieto-Barajas, L. E. and Prünster, I. (2013). Modeling with Normalized Random Measure Mixture Models. Statistical Science. Vol. 28, No. 3, 313-334.
2.- James, L.F., Lijoi, A. and Prünster, I. (2009). Posterior analysis for normalized random measure with independent increments. Scand. J. Statist 36, 76-97.
3.- Kon Kam King, G., Arbel, J. and Prünster, I. (2016). Species Sensitivity Distribution revisited: a Bayesian nonparametric approach. In preparation.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
## Not run: ### Example 1 # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI2cens(x, x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ## Not run: ### Example 2 # Data data(salinity) # Fitting the model under special specifications out <- MixNRMI2cens( xleft = salinity$left, xright = salinity$right, Nit = 5000, distr.pz0 = 10, mu.pz0 = 1, sigma.pz0 = 2 ) # Plotting density estimate + 95% credible interval attach(out) plot(out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) detach() ## End(Not run)
## Not run: ### Example 1 # Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixNRMI2cens(x, x) # Plotting density estimate + 95% credible interval plot(out) ## End(Not run) ## Not run: ### Example 2 # Data data(salinity) # Fitting the model under special specifications out <- MixNRMI2cens( xleft = salinity$left, xright = salinity$right, Nit = 5000, distr.pz0 = 10, mu.pz0 = 1, sigma.pz0 = 2 ) # Plotting density estimate + 95% credible interval attach(out) plot(out) # Plotting number of clusters par(mfrow = c(2, 1)) plot(R, type = "l", main = "Trace of R") hist(R, breaks = min(R - 0.5):max(R + 0.5), probability = TRUE) detach() ## End(Not run)
This function calls the PYdensity function from package BNPmix, to allow fitting a Pitman-Yor process mixture to the data.
MixPY1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Gama = 0.4, asigma = 2, bsigma = 1/var(x), Nx = 100, Nit = 1500, Pbi = 0.5, epsilon = NULL, printtime = TRUE, extras = TRUE )
MixPY1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Gama = 0.4, asigma = 2, bsigma = 1/var(x), Nx = 100, Nit = 1500, Pbi = 0.5, epsilon = NULL, printtime = TRUE, extras = TRUE )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See |
Gama |
Numeric constant. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior on the standard deviation of the mixture kernel. Default value suggested by package BNPmix. |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior on the standard deviation of the mixture kernel. Default value suggested by package BNPmix. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means and weights |
The function returns a MixPY1 object. It is based on a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
R |
Numeric vector of
|
S |
Numeric vector of |
Allocs |
List of |
means |
List of |
weights |
List of
|
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
distr.k |
Integer corresponding to the kernel chosen for the mixture. Always 1, since the Pitman-Yor process is only written to work with Gaussian kernels. |
data |
Data used for the fit |
PY_params |
A named list with the parameters of the Pitman-Yor process |
# Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixPY1(x) # Plotting density estimate + 95% credible interval plot(out)
# Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixPY1(x) # Plotting density estimate + 95% credible interval plot(out)
This function calls the PYdensity function from package BNPmix, to allow fitting a Pitman-Yor process mixture to the data.
MixPY2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Gama = 0.4, asigma = 2, bsigma = 1/var(x), Nx = 100, Nit = 1500, Pbi = 0.5, epsilon = NULL, printtime = TRUE, extras = TRUE )
MixPY2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Gama = 0.4, asigma = 2, bsigma = 1/var(x), Nx = 100, Nit = 1500, Pbi = 0.5, epsilon = NULL, printtime = TRUE, extras = TRUE )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See |
Gama |
Numeric constant. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior on the standard deviation of the mixture kernel. Default value suggested by package BNPmix. |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior on the standard deviation of the mixture kernel. Default value suggested by package BNPmix. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means and weights |
The function returns a MixPY2 object. It is based on a list with the following components:
xx |
Numeric vector. Evaluation grid. |
qx |
Numeric array. Matrix
of dimension |
R |
Numeric vector of
|
Allocs |
List of |
means |
List of |
sigmas |
List of |
weights |
List of
|
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
distr.k |
Integer corresponding to the kernel chosen for the mixture. Always 1, since the Pitman-Yor process is only written to work with Gaussian kernels. |
data |
Data used for the fit |
PY_params |
A named list with the parameters of the Pitman-Yor process |
# Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixPY2(x) # Plotting density estimate + 95% credible interval plot(out)
# Data data(acidity) x <- acidity # Fitting the model under default specifications out <- MixPY2(x) # Plotting density estimate + 95% credible interval plot(out)
Multiple chains of MixNRMI1
multMixNRMI1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
multMixNRMI1( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.p0 |
The distribution name for the centering measure. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure: 1 = Normal; 2 = Gamma; 3 = Beta. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior
on the standard deviation of the mixture kernel |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior
on the standard deviation of the mixture kernel |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling sigma. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
nchains |
The number of chains to run. |
parallel |
Whether to run the chains in parallel. Only works on UNIX-like systems as it rests on Fork parallelism |
ncores |
Number of cores for the parallel run. Defaults to parallel::detectCores(), i.e. the maximum number of cores detected by R on your system. |
a list containing the multiple fits.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
data(acidity) multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2)
data(acidity) multMixNRMI1(acidity, parallel = TRUE, Nit = 10, ncores = 2)
Multiple chains of MixNRMI1cens
multMixNRMI1cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
multMixNRMI1cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.p0 = "normal", asigma = 0.5, bsigma = 0.5, delta_S = 3, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
xleft |
Numeric vector. Lower limit of interval censoring. For exact data the same as xright |
xright |
Numeric vector. Upper limit of interval censoring. For exact data the same as xleft. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.p0 |
The distribution name for the centering measure. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure: 1 = Normal; 2 = Gamma; 3 = Beta. |
asigma |
Numeric positive constant. Shape parameter of the gamma prior
on the standard deviation of the mixture kernel |
bsigma |
Numeric positive constant. Rate parameter of the gamma prior
on the standard deviation of the mixture kernel |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling sigma. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of Nit. |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
nchains |
The number of chains to run. |
parallel |
Whether to run the chains in parallel. Only works on UNIX-like systems as it rests on Fork parallelism |
ncores |
Number of cores for the parallel run. Defaults to parallel::detectCores(), i.e. the maximum number of cores detected by R on your system. |
a list containing the multiple fits.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
data(salinity) multMixNRMI1cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2)
data(salinity) multMixNRMI1cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2)
Multiple chains of MixNRMI2
multMixNRMI2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = FALSE, ncores = parallel::detectCores() )
multMixNRMI2( x, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = FALSE, ncores = parallel::detectCores() )
x |
Numeric vector. Data set to which the density is fitted. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.py0 |
The distribution name for the centering measure for locations. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure for locations: 1 = Normal; 2 = Gamma; 3 = Beta. |
distr.pz0 |
The distribution name for the centering measure for scales. Allowed names are "gamma", or an integer number identifying the centering measure for
scales: 2 = Gamma. For more options use |
mu.pz0 |
Numeric constant. Prior mean of the centering measure for scales. |
sigma.pz0 |
Numeric constant. Prior standard deviation of the centering measure for scales. |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the scales. |
kappa |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the location parameters. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. If 'adaptive=TRUE', 'delta_U'is the starting value for the adaptation. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, sigmas, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
nchains |
The number of chains to run. |
parallel |
Whether to run the chains in parallel. Only works on UNIX-like systems as it rests on Fork parallelism |
ncores |
Number of cores for the parallel run. Defaults to parallel::detectCores(), i.e. the maximum number of cores detected by R on your system. |
a list containing the multiple fits.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
data(acidity) multMixNRMI2(acidity, parallel = TRUE, Nit = 10, ncores = 2)
data(acidity) multMixNRMI2(acidity, parallel = TRUE, Nit = 10, ncores = 2)
Multiple chains of MixNRMI2cens
multMixNRMI2cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
multMixNRMI2cens( xleft, xright, probs = c(0.025, 0.5, 0.975), Alpha = 1, Kappa = 0, Gama = 0.4, distr.k = "normal", distr.py0 = "normal", distr.pz0 = "gamma", mu.pz0 = 3, sigma.pz0 = sqrt(10), delta_S = 4, kappa = 2, delta_U = 2, Meps = 0.01, Nx = 150, Nit = 1500, Pbi = 0.1, epsilon = NULL, printtime = TRUE, extras = TRUE, adaptive = FALSE, nchains = 4, parallel = TRUE, ncores = parallel::detectCores() )
xleft |
Numeric vector. Lower limit of interval censoring. For exact data the same as xright |
xright |
Numeric vector. Upper limit of interval censoring. For exact data the same as xleft. |
probs |
Numeric vector. Desired quantiles of the density estimates. |
Alpha |
Numeric constant. Total mass of the centering measure. See details. |
Kappa |
Numeric positive constant. See details. |
Gama |
Numeric constant. |
distr.k |
The distribution name for the kernel. Allowed names are "normal", "gamma", "beta", "double exponential", "lognormal" or their common abbreviations "norm", "exp", or an integer number identifying the mixture kernel: 1 = Normal; 2 = Gamma; 3 = Beta; 4 = Double Exponential; 5 = Lognormal. |
distr.py0 |
The distribution name for the centering measure for locations. Allowed names are "normal", "gamma", "beta", or their common abbreviations "norm", "exp", or an integer number identifying the centering measure for locations: 1 = Normal; 2 = Gamma; 3 = Beta. |
distr.pz0 |
The distribution name for the centering measure for scales. Allowed names are "gamma", "lognormal", "half-Cauchy", "half-normal", "half-student", "uniform" and "truncated normal", or their common abbreviations "norm", "exp", "lnorm", "halfcauchy", "halfnorm", "halft" and "unif", or an integer number identifying the centering measure for scales: 2 = Gamma, 5 = Lognormal, 6 = Half Cauchy, 7 = Half Normal, 8 = Half Student-t, 9 = Uniform, 10 = Truncated Normal. |
mu.pz0 |
Numeric constant. Prior mean of the centering measure for scales. |
sigma.pz0 |
Numeric constant. Prior standard deviation of the centering measure for scales. |
delta_S |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the scales. |
kappa |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the location parameters. |
delta_U |
Numeric positive constant. Metropolis-Hastings proposal variation coefficient for sampling the latent U. If 'adaptive=TRUE', 'delta_U'is the starting value for the adaptation. |
Meps |
Numeric constant. Relative error of the jump sizes in the continuous component of the process. Smaller values imply larger number of jumps. |
Nx |
Integer constant. Number of grid points for the evaluation of the density estimate. |
Nit |
Integer constant. Number of MCMC iterations. |
Pbi |
Numeric constant. Burn-in period proportion of |
epsilon |
Numeric constant. Extension to the evaluation grid range. See details. |
printtime |
Logical. If TRUE, prints out the execution time. |
extras |
Logical. If TRUE, gives additional objects: means, sigmas, weights and Js. |
adaptive |
Logical. If TRUE, uses an adaptive MCMC strategy to sample the latent U (adaptive delta_U). |
nchains |
The number of chains to run. |
parallel |
Whether to run the chains in parallel. Only works on UNIX-like systems as it rests on Fork parallelism |
ncores |
Number of cores for the parallel run. Defaults to parallel::detectCores(), i.e. the maximum number of cores detected by R on your system. |
a list containing the multiple fits.
MixNRMI2
, MixNRMI1cens
,
MixNRMI2cens
, multMixNRMI1
data(salinity) ## Not run: multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 20, ncores = 2) ## End(Not run)
data(salinity) ## Not run: multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 20, ncores = 2) ## End(Not run)
Determines the jump heights of an increasing additive process by inverting
the M(v) function. Use a truncation level based on expected moments of the NGG process (thresholdGG
).
For internal use.
MvInv(eps, u = 0.5, alpha = 1, kappa = 1, gama = 1/2, N = 3001)
MvInv(eps, u = 0.5, alpha = 1, kappa = 1, gama = 1/2, N = 3001)
eps |
Dummy argument kept for consistency with past versions of the functions |
u |
Real number. The value of the latent variable at the current step. |
alpha |
Numeric constant. Total mass of the centering measure. |
kappa |
Numeric positive constant. |
gama |
Numeric constant. Discount parameter of the NRMI process. |
N |
Number of steps in the discretization scheme for the grid inversion. ## The function has been optimised but it is morally defined as: function(eps, u = 0.5, alpha = 1, kappa = 1, gama = 1 / 2, N = 3001) n <- length(w) v <- rep(NA, n) x <- -log(seq(from = exp(-1e-05), to = exp(-10), length = N)) f <- alpha / gamma(1 - gama) * x^(-(1 + gama)) * exp(-(u + kappa) * x) dx <- diff(x) h <- (f[-1] + f[-N]) / 2 Mv <- rep(0, N) for (i in seq(N - 1, 1)) Mv[i] <- Mv[i + 1] + dx[i] * h[i] for (j in seq(n)) v[j] <- x[which.min(Mv > w[j])] return(v) |
This is a function to visualize the clustering induced by the BNP model. The data points are plotted with a color reflecting their cluster.
plot_clustering_and_CDF(fit, clustering, label_vector = NULL)
plot_clustering_and_CDF(fit, clustering, label_vector = NULL)
fit |
The fitted object, obtained from one of the MixNRMIx functions |
clustering |
A vector of integers with the same length as the data, representing the allocation variable for data each point. |
label_vector |
A vector of data labels to be plotted, to provide some identification to each point. |
A plot of the Cumulative Distribution Function (or Turnbull estimate for censored data) with data points whose color denotes the cluster allocation. For censored data, right or left censored data points are not represented, while interval censored data points are represented at the middle of the censoring interval.
This plots the prior distribution on the number of components for the stable process. The Dirichlet process is provided for comparison.
plot_prior_number_of_components( n, Gama, Alpha = 1, grid = NULL, silence = TRUE )
plot_prior_number_of_components( n, Gama, Alpha = 1, grid = NULL, silence = TRUE )
n |
Number of data points |
Gama |
Numeric constant. 0 <= Gama <=1. |
Alpha |
Numeric constant. Total mass of the centering measure for the Dirichlet process. |
grid |
Integer vector. Level of truncation when computing the expectation. Defaults to n. If greater than n, it is fixed to n. |
silence |
Boolean. Whether to print the current calculation step for the Stable process, as the function can be long |
A plot with the prior distribution on the number of components.
plot_prior_number_of_components(50, 0.4)
plot_prior_number_of_components(50, 0.4)
The density estimate is the mean posterior density computed on the data points.
## S3 method for class 'multNRMI' plot(x, ...)
## S3 method for class 'multNRMI' plot(x, ...)
x |
An object of class multNRMI |
... |
Further arguments to be passed to generic functions, ignored at the moment |
A graph with the density estimate, the 95% credible interval. Includes a histogram if the data is non censored.
data(salinity) fit <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) plot(fit)
data(salinity) fit <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) plot(fit)
The density estimate is the mean posterior density computed on the data points.
## S3 method for class 'NRMI1' plot(x, ...)
## S3 method for class 'NRMI1' plot(x, ...)
x |
A fitted object of class NRMI1 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A graph with the density estimate, the 95% credible interval and a histogram of the data
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out) ## Example for censored data data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, Nit = 50) plot(out)
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out) ## Example for censored data data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, Nit = 50) plot(out)
The density estimate is the mean posterior density computed on the data points.
## S3 method for class 'NRMI2' plot(x, ...)
## S3 method for class 'NRMI2' plot(x, ...)
x |
A fitted object of class NRMI2 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A graph with the density estimate, the 95% credible interval and a histogram of the data
## Example for non censored data data(acidity) out <- MixNRMI2(acidity, Nit = 20) plot(out) ## Example for censored data data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) plot(out)
## Example for non censored data data(acidity) out <- MixNRMI2(acidity, Nit = 20) plot(out) ## Example for censored data data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) plot(out)
Plot the density estimate and the 95% credible interval
## S3 method for class 'PY1' plot(x, ...)
## S3 method for class 'PY1' plot(x, ...)
x |
A fitted object of class PY1 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A graph with the density estimate, the 95% credible interval and a histogram of the data
data(acidity) out <- MixPY1(acidity, Nit = 50) plot(out)
data(acidity) out <- MixPY1(acidity, Nit = 50) plot(out)
Plot the density estimate and the 95% credible interval
## S3 method for class 'PY2' plot(x, ...)
## S3 method for class 'PY2' plot(x, ...)
x |
A fitted object of class PY2 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A graph with the density estimate, the 95% credible interval and a histogram of the data
data(acidity) out <- MixPY2(acidity, Nit = 50) plot(out)
data(acidity) out <- MixPY2(acidity, Nit = 50) plot(out)
Plot the Turnbull CDF and fitted CDF for censored data.
plotCDF_censored(fit)
plotCDF_censored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1cens or MixNRMI2cens. |
Plot of the empirical and fitted CDF for non censored data.
set.seed(150520) data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) BNPdensity:::plotCDF_censored(out)
set.seed(150520) data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, extras = TRUE, Nit = 100) BNPdensity:::plotCDF_censored(out)
Plot the empirical and fitted CDF for non censored data.
plotCDF_noncensored(fit)
plotCDF_noncensored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2. |
Plot of the empirical and fitted CDF for non censored data.
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 10) BNPdensity:::plotCDF_noncensored(out)
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 10) BNPdensity:::plotCDF_noncensored(out)
The density estimate is the mean posterior density computed on the data points. It is not possible to display a histogram for censored data.
plotfit_censored(fit)
plotfit_censored(fit)
fit |
A fitted object of class NRMI1cens or NRMI2cens |
A graph with the density estimate and the 95% credible interval
data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out)
data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out)
The density estimate is the mean posterior density computed on the data points.
plotfit_noncensored(fit)
plotfit_noncensored(fit)
fit |
A fitted object of class NRMI1 or NRMI2 |
A graph with the density estimate, the 95% credible interval and a histogram of the data
data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out)
data(acidity) out <- MixNRMI1(acidity, Nit = 50) plot(out)
Plot the density for censored data.
plotPDF_censored(fit)
plotPDF_censored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1cens or MixNRMI2cens. |
Plot of the density and a histogram for non censored data.
set.seed(150520) data(salinity) out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, extras = TRUE, Nit = 100) BNPdensity:::plotPDF_censored(out)
set.seed(150520) data(salinity) out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, extras = TRUE, Nit = 100) BNPdensity:::plotPDF_censored(out)
Plot the density and a histogram for non censored data.
plotPDF_noncensored(fit)
plotPDF_noncensored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2. |
Plot of the density and a histogram for non censored data.
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::plotPDF_noncensored(out)
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::plotPDF_noncensored(out)
Plot the percentile-percentile graph for non censored data, using the Turnbull estimator the position of the percentiles.
pp_plot_censored(fit)
pp_plot_censored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1cens or MixNRMI2cens. |
Percentile-percentile graph using the Turnbull estimator
set.seed(150520) data(salinity) out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, extras = TRUE, Nit = 100) BNPdensity:::pp_plot_censored(out)
set.seed(150520) data(salinity) out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, extras = TRUE, Nit = 100) BNPdensity:::pp_plot_censored(out)
Plot the percentile-percentile graph for non censored data.
pp_plot_noncensored(fit)
pp_plot_noncensored(fit)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2. |
Percentile-percentile plot for non censored data.
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::pp_plot_noncensored(out)
set.seed(150520) data(acidity) out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) BNPdensity:::pp_plot_noncensored(out)
S3 method for class 'multNRMI'
## S3 method for class 'multNRMI' print(x, ...)
## S3 method for class 'multNRMI' print(x, ...)
x |
An object of class multNRMI |
... |
Further arguments to be passed to generic functions, ignored at the moment |
A visualization of the important information about the object
data(salinity) out <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) print(out)
data(salinity) out <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) print(out)
S3 method for class 'MixNRMI1'
## S3 method for class 'NRMI1' print(x, ...)
## S3 method for class 'NRMI1' print(x, ...)
x |
A fitted object of class NRMI1 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A visualization of the important information about the object
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) print(out) ## Example for censored data data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, Nit = 50) print(out)
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) print(out) ## Example for censored data data(salinity) out <- MixNRMI1cens(salinity$left, salinity$right, Nit = 50) print(out)
S3 method for class 'MixNRMI2'
## S3 method for class 'NRMI2' print(x, ...)
## S3 method for class 'NRMI2' print(x, ...)
x |
A fitted object of class NRMI2 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A visualization of the important information about the object
#' ## Example for censored data data(acidity) out <- MixNRMI2(acidity, Nit = 20) print(out) data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) print(out)
#' ## Example for censored data data(acidity) out <- MixNRMI2(acidity, Nit = 20) print(out) data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) print(out)
S3 method for class 'PY1'
## S3 method for class 'PY1' print(x, ...)
## S3 method for class 'PY1' print(x, ...)
x |
A fitted object of class PY1 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A visualization of the important information about the object
## Example for non censored data data(acidity) out <- MixPY1(acidity, Nit = 50) print(out)
## Example for non censored data data(acidity) out <- MixPY1(acidity, Nit = 50) print(out)
S3 method for class 'PY2'
## S3 method for class 'PY2' print(x, ...)
## S3 method for class 'PY2' print(x, ...)
x |
A fitted object of class PY2 |
... |
Further arguments to be passed to generic function, ignored at the moment |
A visualization of the important information about the object
## Example for non censored data data(acidity) out <- MixPY2(acidity, Nit = 50) print(out)
## Example for non censored data data(acidity) out <- MixPY2(acidity, Nit = 50) print(out)
This function is intended to help with compatibility with the previous versions of the package.
process_dist_name(distname)
process_dist_name(distname)
distname |
Can be an integer or a distribution name. Allowed names are "normal", "gamma", "beta", "exponential", "lognormal", "half-Cauchy", "half-normal", "half-student", "uniform" and "truncated normal", or their common abbreviations "norm", "exp", "halfcauchy", "halfnorm", "halft" and "unif". |
an integer both if distname is an integer or a character
This function may be rather slow for many iterations/many data because it relies on numerical inversion of the mixture Cumulative Distribution Function. set.seed(150520) data(salinity) out <- MixNRMI1cens(xleft = salinity$left, xright = salinity$right, extras = TRUE, Nit = 100) BNPdensity:::qq_plot_censored(out)
qq_plot_censored(fit, thinning_to = 500)
qq_plot_censored(fit, thinning_to = 500)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2, MixMRMI1cens or MixMRMI2cens |
thinning_to |
How many iterations to compute the mean posterior quantiles |
quantile-quantile plot for non censored data.
This function may be rather slow for many iterations/many data because it relies on numerical inversion of the mixture Cumulative Distribution Function.
qq_plot_noncensored(fit, thinning_to = 500)
qq_plot_noncensored(fit, thinning_to = 500)
fit |
The result of the fit, obtained through the function MixNRMI1 or MixNRMI2, MixMRMI1cens or MixMRMI2cens |
thinning_to |
How many iterations to compute the mean posterior quantiles |
quantile-quantile plot for non censored data.
### Not run # set.seed(150520) # data(acidity) # out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) # BNPdensity:::qq_plot_noncensored(out)
### Not run # set.seed(150520) # data(acidity) # out <- MixNRMI1(acidity, extras = TRUE, Nit = 100) # BNPdensity:::qq_plot_noncensored(out)
72-hour acute salinity tolerance (LC50 values) of riverine macro-invertebrates.
A data frame with 108 observations on the following two variables:
A numeric vector.
A numeric vector.
fitdistrplus
R-package
Kefford, B.J., Nugegoda, D., Metzeling, L., Fields, E. 2006. Validating species sensitivity distributions using salinity tolerance of riverine macroinvertebrates in the southern Murray-darling Basin (Victoria, Australia). Canadian Journal of Fisheries and Aquatic Science, 63, 1865-1877.
data(salinity) hist(salinity$left)
data(salinity) hist(salinity$left)
S3 method for class 'multNRMI'
## S3 method for class 'multNRMI' summary(object, number_of_clusters = FALSE, ...)
## S3 method for class 'multNRMI' summary(object, number_of_clusters = FALSE, ...)
object |
A fitted object of class NRMI1cens |
number_of_clusters |
Whether to compute the optimal number of clusters, which can be a time-consuming operation (see |
... |
Further arguments to be passed to generic function, ignored at the moment |
Prints out the text for the summary S3 methods
data(salinity) out <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) summary(out)
data(salinity) out <- multMixNRMI2cens(salinity$left, salinity$right, parallel = TRUE, Nit = 10, ncores = 2) summary(out)
S3 method for class 'MixNRMI1'
## S3 method for class 'NRMI1' summary(object, number_of_clusters = FALSE, ...)
## S3 method for class 'NRMI1' summary(object, number_of_clusters = FALSE, ...)
object |
A fitted object of class NRMI1 |
number_of_clusters |
Whether to compute the optimal number of clusters, which can be a time-consuming operation (see |
... |
Further arguments to be passed to generic function, ignored at the moment |
Prints out the text for the summary S3 methods
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) summary(out)
## Example for non censored data data(acidity) out <- MixNRMI1(acidity, Nit = 50) summary(out)
S3 method for class 'MixNRMI2'
## S3 method for class 'NRMI2' summary(object, number_of_clusters = FALSE, ...)
## S3 method for class 'NRMI2' summary(object, number_of_clusters = FALSE, ...)
object |
A fitted object of class NRMI2 |
number_of_clusters |
Whether to compute the optimal number of clusters, which can be a time-consuming operation (see |
... |
Further arguments to be passed to generic function, ignored at the moment |
Prints out the text for the summary S3 methods
data(acidity) out <- MixNRMI2(acidity, Nit = 20) summary(out) data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) summary(out)
data(acidity) out <- MixNRMI2(acidity, Nit = 20) summary(out) data(salinity) out <- MixNRMI2cens(salinity$left, salinity$right, Nit = 20) summary(out)
S3 method for class 'PY1'
## S3 method for class 'PY1' summary(object, number_of_clusters = FALSE, ...)
## S3 method for class 'PY1' summary(object, number_of_clusters = FALSE, ...)
object |
A fitted object of class PY1 |
number_of_clusters |
Whether to compute the optimal number of clusters, which can be a time-consuming operation (see |
... |
Further arguments to be passed to generic function, ignored at the moment |
Prints out the text for the summary S3 methods
## Example for non censored data data(acidity) out <- MixPY1(acidity, Nit = 50) summary(out)
## Example for non censored data data(acidity) out <- MixPY1(acidity, Nit = 50) summary(out)
S3 method for class 'PY2'
## S3 method for class 'PY2' summary(object, number_of_clusters = FALSE, ...)
## S3 method for class 'PY2' summary(object, number_of_clusters = FALSE, ...)
object |
A fitted object of class PY2 |
number_of_clusters |
Whether to compute the optimal number of clusters, which can be a time-consuming operation (see |
... |
Further arguments to be passed to generic function, ignored at the moment |
Prints out the text for the summary S3 methods
## Example for non censored data data(acidity) out <- MixPY2(acidity, Nit = 50) summary(out)
## Example for non censored data data(acidity) out <- MixPY2(acidity, Nit = 50) summary(out)
Common text for the summary S3 methods
summarytext( fit, kernel_comment, BNP_process_comment, number_of_clusters = FALSE )
summarytext( fit, kernel_comment, BNP_process_comment, number_of_clusters = FALSE )
fit |
NRMIx or PYx object |
kernel_comment |
Text specific to the parametric and nonparametric nature of the model |
BNP_process_comment |
Text specific to the nonparametric process, NRMI or Pitman-Yor |
number_of_clusters |
Flag to decide whether to compute the optimal clustering |
Prints out the text for the summary S3 methods
This is a convenience function which works when coda is not yet loaded by the user. If coda is loaded, it gets masked. See also file multMixNRMI.R
traceplot(fitlist)
traceplot(fitlist)
fitlist |
Output of multMixNRMI. |
A traceplot for multiple chains.