rnegbinRw {bayesm} | R Documentation

## MCMC Algorithm for Negative Binomial Regression

### Description

`rnegbinRw` implements a Random Walk Metropolis Algorithm for the Negative Binomial (NBD) regression model. beta | alpha and alpha | beta are drawn with two different random walks.

### Usage

    rnegbinRw(Data, Prior, Mcmc)

### Arguments

`Data` | list(y,X)
`Prior` | list(betabar,A,a,b)
`Mcmc` | list(R,keep,s_beta,s_alpha,beta0

### Details

Model: y ~ NBD(mean=lambda, over-dispersion=alpha).
lambda=exp(x'beta)

Prior: beta ~ N(betabar,A^{-1})
alpha ~ Gamma(a,b).
note: prior mean of alpha = a/b, variance = a/(b^2)

list arguments contain:

`y`
  nobs vector of counts (0,1,2,...)
`X`
  nobs x nvar matrix
`betabar`
  nvar x 1 prior mean (def: 0)
`A`
  nvar x nvar pds prior prec matrix (def: .01I)
`a`
  Gamma prior parm (def: .5)
`b`
  Gamma prior parm (def: .1)
`R`
  number of MCMC draws
`keep`
  MCMC thinning parm: keep every keepth draw (def: 1)
`s\_beta`
  scaling for beta| alpha RW inc cov matrix (def: 2.93/sqrt(nvar)
`s\_alpha`
  scaling for alpha | beta RW inc cov matrix (def: 2.93)

### Value

a list containing:

`betadraw` | R/keep x nvar array of beta draws
`alphadraw` | R/keep vector of alpha draws
`llike` | R/keep vector of log-likelihood values evaluated at each draw
`acceptrbeta` | acceptance rate of the beta draws
`acceptralpha` | acceptance rate of the alpha draws

### Note

The NBD regression encompasses Poisson regression in the sense that as alpha goes to infinity the NBD distribution tends toward the Poisson.
For "small" values of alpha, the dependent variable can be extremely variable so that a large number of observations may be required to obtain precise inferences.

### Author(s)

Sridhar Narayanam & Peter Rossi, Graduate School of Business, University of Chicago, Peter.Rossi@ChicagoGsb.edu.

### References

For further discussion, see Bayesian Statistics and Marketing by Rossi, Allenby, McCulloch.
http://gsbwww.uchicago.edu/fac/peter.rossi/research/bsm.html

### See Also

`rhierNegbinRw`

### Examples

    ##
    if(nchar(Sys.getenv("LONG_TEST")) != 0)  {R=1000} else {R=10}

    set.seed(66)
    simnegbin =
    function(X, beta, alpha) {
    #   Simulate from the Negative Binomial Regression
    lambda = exp(X %*% beta)
    y=NULL
    for (j in 1:length(lambda))
        y = c(y,rnbinom(1,mu = lambda[j],size = alpha))
    return(y)
    }

    nobs = 500
    nvar=2            # Number of X variables
    alpha = 5
    Vbeta = diag(nvar)*0.01

    # Construct the regdata (containing X)
    simnegbindata = NULL
    beta = c(0.6,0.2)
    X = cbind(rep(1,nobs),rnorm(nobs,mean=2,sd=0.5))
    simnegbindata = list(y=simnegbin(X,beta,alpha), X=X, beta=beta)
    Data = simnegbindata
    betabar = rep(0,nvar)
    A = 0.01 * diag(nvar)
    a = 0.5; b = 0.1
    Prior = list(betabar=betabar, A=A, a=a, b=b)

    keep =1
    s_beta=2.93/sqrt(nvar); s_alpha=2.93
    Mcmc = list(R=R, keep = keep, s_beta=s_beta, s_alpha=s_alpha)
    out = rnegbinRw(Data, Prior, Mcmc)

    cat(" betadraws ",fill=TRUE)
    mat=apply(out$betadraw,2,quantile,probs=c(.01,.05,.5,.95,.99))
    mat=rbind(beta,mat); rownames(mat)[1]="beta"; print(mat)
    cat(" alphadraws ",fill=TRUE)
    mat=apply(matrix(out$alphadraw),2,quantile,probs=c(.01,.05,.5,.95,.99))
    mat=rbind(as.vector(alpha),mat); rownames(mat)[1]="alpha"; print(mat)

[Package bayesm version 2.0-8 Index]