rivDP {bayesm} R Documentation

## Linear "IV" Model with DP Process Prior for Errors

### Description

rivDP is a Gibbs Sampler for a linear structural equation with an arbitrary number of instruments. rivDP uses a mixture-of-normals for the structural and reduced form equations implemented with a Dirichlet Process prior.

### Usage

rivDP(Data, Prior, Mcmc)

### Arguments

 Data list(y, x, w, z) Prior list(md, Ad, mbg, Abg, lambda, Prioralpha, lambda_hyper) Mcmc list(R, keep, nprint, maxuniq, SCALE, gridsize)

### Details

#### Model and Priors

x = z'\delta + e1
y = \beta*x + w'\gamma + e2
e1,e2 \sim N(\theta_{i}) where \theta_{i} represents \mu_{i}, \Sigma_{i}

Note: Error terms have non-zero means. DO NOT include intercepts in the z or w matrices. This is different from rivGibbs which requires intercepts to be included explicitly.

\delta \sim N(md, Ad^{-1})
vec(\beta, \gamma) \sim N(mbg, Abg^{-1})
\theta_{i} \sim G
G \sim DP(alpha, G_0)

alpha \sim (1-(alpha-alpha_{min})/(alpha_{max}-alpha{min}))^{power}
where alpha_{min} and alpha_{max} are set using the arguments in the reference below. It is highly recommended that you use the default values for the hyperparameters of the prior on alpha.

G_0 is the natural conjugate prior for (\mu,\Sigma): \Sigma \sim IW(nu, vI) and \mu|\Sigma \sim N(0, \Sigma(x) a^{-1})
These parameters are collected together in the list \lambda. It is highly recommended that you use the default settings for these hyper-parameters.

\lambda(a, nu, v):
a \sim uniform[alim[1], alimb[2]]
nu \sim dim(data)-1 + exp(z)
z \sim uniform[dim(data)-1+nulim[1], nulim[2]]
v \sim uniform[vlim[1], vlim[2]]

#### Argument Details

Data = list(y, x, w, z)

 y:  n x 1 vector of obs on LHS variable in structural equation x:  n x 1 vector of obs on "endogenous" variable in structural equation w:  n x j matrix of obs on "exogenous" variables in the structural equation z:  n x p matrix of obs on instruments

Prior = list(md, Ad, mbg, Abg, lambda, Prioralpha, lambda_hyper) [optional]

 md:  p-length prior mean of delta (def: 0) Ad:  p x p PDS prior precision matrix for prior on delta (def: 0.01*I) mbg:  (j+1)-length prior mean vector for prior on beta,gamma (def: 0) Abg:  (j+1)x(j+1) PDS prior precision matrix for prior on beta,gamma (def: 0.01*I) Prioralpha: list(Istarmin, Istarmax, power) $Istarmin:  is expected number of components at lower bound of support of alpha (def: 1) $Istarmax:  is expected number of components at upper bound of support of alpha (def: floor(0.1*length(y))) $power:  is the power parameter for alpha prior (def: 0.8) lambda_hyper: list(alim, nulim, vlim) $alim:  defines support of a distribution (def: c(0.01, 10)) $nulim:  defines support of nu distribution (def: c(0.01, 3)) $vlim:  defines support of v distribution (def: c(0.1, 4))

Mcmc = list(R, keep, nprint, maxuniq, SCALE, gridsize) [only R required]

 R:  number of MCMC draws keep:  MCMC thinning parameter: keep every keepth draw (def: 1) nprint:  print the estimated time remaining for every nprint'th draw (def: 100, set to 0 for no print) maxuniq:  storage constraint on the number of unique components (def: 200) SCALE:  scale data (def: TRUE) gridsize:  gridsize parameter for alpha draws (def: 20)

#### nmix Details

nmix is a list with 3 components. Several functions in the bayesm package that involve a Dirichlet Process or mixture-of-normals return nmix. Across these functions, a common structure is used for nmix in order to utilize generic summary and plotting functions.

 probdraw: ncomp x R/keep matrix that reports the probability that each draw came from a particular component (here, a one-column matrix of 1s) zdraw:  R/keep x nobs matrix that indicates which component each draw is assigned to (here, null) compdraw: A list of R/keep lists of ncomp lists. Each of the inner-most lists has 2 elemens: a vector of draws for mu and a matrix of draws for the Cholesky root of Sigma.

### Value

A list containing:

 deltadraw  R/keep x p array of delta draws betadraw  R/keep x 1 vector of beta draws alphadraw  R/keep x 1 vector of draws of Dirichlet Process tightness parameter Istardraw  R/keep x 1 vector of draws of the number of unique normal components gammadraw  R/keep x j array of gamma draws nmix  a list containing: probdraw, zdraw, compdraw (see “nmix Details” section)

### Author(s)

Peter Rossi, Anderson School, UCLA, perossichi@gmail.com.

### References

For further discussion, see "A Semi-Parametric Bayesian Approach to the Instrumental Variable Problem," by Conley, Hansen, McCulloch and Rossi, Journal of Econometrics (2008).

See also, Chapter 4, Bayesian Non- and Semi-parametric Methods and Applications by Peter Rossi.

rivGibbs

### Examples

if(nchar(Sys.getenv("LONG_TEST")) != 0) {R=2000} else {R=10}
set.seed(66)

## simulate scaled log-normal errors and run
k = 10
delta = 1.5
Sigma = matrix(c(1, 0.6, 0.6, 1), ncol=2)
N = 1000
tbeta = 4
scalefactor = 0.6
root = chol(scalefactor*Sigma)
mu = c(1,1)

## compute interquartile ranges
ninterq = qnorm(0.75) - qnorm(0.25)
error = matrix(rnorm(100000*2), ncol=2)%*%root
error = t(t(error)+mu)
Err = t(t(exp(error))-exp(mu+0.5*scalefactor*diag(Sigma)))
lnNinterq = quantile(Err[,1], prob=0.75) - quantile(Err[,1], prob=0.25)

## simulate data
error = matrix(rnorm(N*2), ncol=2)%*%root
error = t(t(error)+mu)
Err = t(t(exp(error))-exp(mu+0.5*scalefactor*diag(Sigma)))

## scale appropriately
Err[,1] = Err[,1]*ninterq/lnNinterq
Err[,2] = Err[,2]*ninterq/lnNinterq
z = matrix(runif(k*N), ncol=k)
x = z%*%(delta*c(rep(1,k))) + Err[,1]
y = x*tbeta + Err[,2]

## specify data input and mcmc parameters
Data = list();
Data$z = z Data$x = x
Data$y = y Mcmc = list() Mcmc$maxuniq = 100
Mcmc$R = R end = Mcmc$R

out = rivDP(Data=Data, Mcmc=Mcmc)

cat("Summary of Beta draws", fill=TRUE)
summary(out$betadraw, tvalues=tbeta) ## plotting examples if(0){ plot(out$betadraw, tvalues=tbeta)
plot(out\$nmix)  # plot "fitted" density of the errors
}


[Package bayesm version 3.1-6 Index]