hurdlr (version 0.1)

hurdle: Hurdle Model Count Data Regression

Description

hurdle is used to fit single or double-hurdle regression models to count data via Bayesian inference.

Usage

hurdle(y, x = NULL, hurd = Inf, dist = c("poisson", "nb", "lognormal"),
  dist.2 = c("gpd", "poisson", "lognormal", "nb"),
  control = hurdle_control(), iters = 1000, burn = 500, nthin = 1,
  plots = FALSE, progress.bar = TRUE)

Arguments

y

numeric response vector.

x

numeric predictor matrix.

hurd

numeric threshold for 'extreme' observations of two-hurdle models. Inf for one-hurdle models.

dist

character specification of response distribution.

dist.2

character specification of response distribution for 'extreme' observations of two-hurdle models.

control

list of parameters for controlling the fitting process, specified by hurdle_control.

iters

number of iterations for the Markov chain to run.

burn

numeric burn-in length.

nthin

numeric thinning rate.

plots

logical operator. TRUE to output plots.

progress.bar

logical operator. TRUE to print progress bar.

Value

hurdle returns a list which includes the items

pD

measure of model dimensionality \(p_D\) where \(p_D = \bar{D} - D(\bar{\theta}\)) is the \("mean posterior deviance - deviance of posterior means"\)

DIC

Deviance Information Criterion where \(DIC = \bar{D} - p_D\)

PPO

Posterior Predictive Ordinate (PPO) measure of fit

CPO

Conditional Predictive Ordinate (CPO) measure of fit

pars.means

posterior mean(s) of third-component parameter(s) if hurd != Inf

ll.means

posterior means of the log-likelihood distributions of all model components

beta.means

posterior means regression coefficients

dev

posterior deviation where \(D = -2LogL\)

beta

posterior distributions of regression coefficients

pars

posterior distribution(s) of third-component parameter(s) if hurd != Inf

Details

Setting dist and dist.2 to be the same distribution creates a single dist-hurdle model, not a double-hurdle model. However, this is being considered in future package updates.

Examples

#Generate some data:
p=0.5; q=0.25; lam=3;
mu=10; sigma=7; xi=0.75;
n=200

set.seed(2016)
y <- rbinom(n,1,p)
nz <- sum(1-y)
extremes <- rbinom(sum(y),1,q)
ne <- sum(extremes)
nt <- n-nz-ne
yt <- sample(mu-1,nt,replace=TRUE,prob=dpois(1:(mu-1),3)/(ppois(mu-1,lam)-ppois(0,lam)))
yz <- round(rgpd(nz,mu,sigma,xi))
y[y==1] <- c(yt,yz)
g <- hurdle(y)