dmixnegbinom: The Negative Binomial Mixture Distribution

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/mederrRank.R

Description

Density function for a mixture of two negative binomial distributions.

Usage

1
dmixnegbinom(x, theta, E, log.p = FALSE)

Arguments

x

vector of (non-negative integer) quantiles.

theta

vector of parameters for the negative binomial distribution mixture.

E

vector of (non-negative integer) expected counts.

log.p

logical; if TRUE, probabilities p are given as log(p).

Details

The mixture of two negative binomial distributions has density

P(N = x) = theta[5] f(x;theta[1],theta[2],E) + (1 - theta[5]) f(x;theta[3],theta[4],E),

where

f(x;α,β,E) = \frac{Γ(α + x)}{Γ(α) x!}\frac{1}{(1 + β/E)^x}\frac{1}{(1 + E/β)^α}

for x = 0,1,…,α, α, β, E > 0 and 0 < theta[5] ≤q 1. The mixture of two negative binomial distributions represents the marginal distribution of the counts N coming from Poisson data with parameter λ and a mixture of two gamma distributions as its prior. For details see the paper by Dumouchel (1999).

Value

dmixnegbinom gives the density corresponding to the E and theta values provided.

Author(s)

Sergio Venturini [email protected],

Jessica A. Myers [email protected]

References

DuMouchel W. (1999), "Bayesian Data Mining in Large Frequency Tables, with an Application to the FDA Spontaneous Reporting System". The American Statistician, 53, 177-190.

Myers, J. A., Venturini, S., Dominici, F. and Morlock, L. (2011), "Random Effects Models for Identifying the Most Harmful Medication Errors in a Large, Voluntary Reporting Database". Technical Report.

See Also

dnbinom, rmixnegbinom.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
## Not run: 
data("simdata", package = "mederrRank")
ni <- simdata@numi
theta0 <- c(10, 6, 100, 100, .1)
ans <- mixnegbinom.em(simdata, theta0, 50000, 0.01,
	se = FALSE, stratified = TRUE)
theta <- ans$theta.hat
N.E <- cbind(ans$N[1:ni], ans$E[1:ni])[sort(ans$N[1:ni], index.return = TRUE)$ix, ]
N.ix <- match(unique(N.E[, 1]), N.E[, 1])
N <- N.E[N.ix, 1]
E <- N.E[N.ix, 2]
dens <- dmixnegbinom(N, theta, E)
hist(N.E[, 1], breaks = 40, freq = FALSE)
points(N, dens)

## End(Not run)

mederrRank documentation built on May 30, 2017, 2:55 a.m.