rnbinom.gf: Generate Time Series with Negative Binomial Distribution and...
In spass: Study Planning and Adaptation of Sample Size

Description Usage Arguments Details Value Source References Examples

rnbinom.gf generates one or more independent time series following the Gamma frailty model. The generated data has negative binomial marginal distribution and the underlying multivariate Gamma frailty an autoregressive covariance structure.

1	rnbinom.gf(n, size, mu, rho, tp)

`n`	number of observations. If length(n) > 1, the length is taken to be the number required.
`size`	dispersion parameter (the shape parameter of the gamma mixing distribution). Must be strictly positive, need not be integer.
`mu`	vector of means of time points: see 'Details'.
`rho`	correlation coefficient of the underlying autoregressive Gamma frailty. Must be between 0 and 1.
`tp`	number of observed time points.

The generated marginal negative binomial distribution with mean mu = μ and size = η has density

(μ/(μ+η))^x Γ(x + η)/(Γ(x+1)Γ(η)) (η/(μ+η))^η

for 0 < μ, 0 < η and x=0, 1, 2, .... Hereby, each entry of vector mu corresponds to one time point. Therefore, each timepoint can have its distinct mean.

Within the Gamma frailty model, the correlation between two frailties of time points t and s for rho = ρ is given by

ρ^|t-s|

for 0 ≤ ρ ≤ 1. Note: this does not correspond to the correlation of observations.

rnbinom.gf returns a matrix of dimension n x tp with marginal negative binomial distribution with means mu, common dispersion parameter size and a correlation induce by the autoregressive multivariate Gamma frailty.

rnbinom.gf computes observations from a Gamma frailty model by Fiocco et. al. 2009 using code contributed by Thomas Asendorf.

Fiocco M, Putter H, Van Houwelingen JC, (2009), A new serially correlated gamma-frailty process for longitudinal count data Biostatistics Vol. 10, No. 2, pp. 245-257.

set.seed(8)
random<-rnbinom.gf(n=1000, size=0.6, mu=1:6, rho=0.8, tp=6)
cor(random)

#Check the marginal distribution of time point 3
plot(table(random[,3])/1000, xlab="Probability", ylab="Observation")
lines(0:26, dnbinom(0:26, mu=3, size=0.6), col="red")
legend("topright",legend=c("Theoretical Marginal Distribution", "Observed Distribution"),
  col=c("red", "black"), lty=1, lwd=c(1,2))