sandwich2 {spass} | R Documentation |
Generate Time Series with Negative Binomial Distribution and Autoregressive Correlation Structure of Order One: NB-INAR(1)
Description
rnbinom.inar1
generates one or more independent time series following the NB-INAR(1) model. The generated data has negative binomial marginal distribution and an autoregressive covariance structure.
Usage
sandwich2(sigma, rho, theta, k, Time, dropout, Model)
Arguments
sigma |
assymptotic standard deviation for Full and subpupulation |
rho |
correlation coefficient of the underlying autoregressive correlation structure. Must be between 0 and 1. |
theta |
correlation absorption coefficient if tinepoints are farther appart |
k |
sample size allocation factor between groups: see 'Details'. |
Time |
vector of measured timepoints |
dropout |
vector describing the percentage of dropout in every timepoint |
Model |
either 1 or 2, describing if 4-regressor or 3-regressor model was used. |
Details
The generated marginal negative binomial distribution with mean mu
= \mu
and size
= \eta
has density
(\mu/(\mu+\eta))^x \Gamma(x + \eta)/(\Gamma(x+1)\Gamma(\eta)) (\eta/(\mu+\eta))^\eta
for 0 < \mu
, 0 < \eta
and x=0, 1, 2, ...
.
Within the NB-INAR(1) model, the correlation between two time points t
and s
for rho
= \rho
is given through
\rho^|t-s|
for 0 \le \rho \le 1
.
Value
rnbinom.inar1
returns a matrix of dimension n
x tp
with marginal negative binomial
distribution with mean mu
and dispersion parameter size
, and an autoregressive correlation structure
between time points.
Source
rnbinom.inar1
computes a reparametrization of the NB-INAR(1) model by McKenzie 1986 using code contributed by Thomas Asendorf.
References
McKenzie Ed (1986), Autoregressive Moving-Average Processes with Negative-Binomial and Geometric Marginal Distributions. Advances in Applied Probability Vol. 18, No. 3, pp. 679-705.
Examples
set.seed(8)
random<-rnbinom.inar1(n=1000, size=0.6, mu=2, 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=2, size=0.6), col="red")
legend("topright",legend=c("Theoretical Marginal Distribution", "Observed Distribution"),
col=c("red", "black"), lty=1, lwd=c(1,2))