Learn R Programming

bayesm (version 3.1-6)

rhierLinearModel: Gibbs Sampler for Hierarchical Linear Model with Normal Heterogeneity

Description

rhierLinearModel implements a Gibbs Sampler for hierarchical linear models with a normal prior.

Usage

rhierLinearModel(Data, Prior, Mcmc)

Value

A list containing:

betadraw

\(nreg x nvar x R/keep\) array of individual regression coef draws

taudraw

\(R/keep x nreg\) matrix of error variance draws

Deltadraw

\(R/keep x nz*nvar\) matrix of Deltadraws

Vbetadraw

\(R/keep x nvar*nvar\) matrix of Vbeta draws

Arguments

Data

list(regdata, Z)

Prior

list(Deltabar, A, nu.e, ssq, nu, V)

Mcmc

list(R, keep, nprint)

Author

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

Details

Model and Priors

nreg regression equations with \(nvar\) \(X\) variables in each equation
\(y_i = X_i\beta_i + e_i\) with \(e_i\) \(\sim\) \(N(0, \tau_i)\)

\(\tau_i\) \(\sim\) nu.e*\(ssq_i/\chi^2_{nu.e}\) where \(\tau_i\) is the variance of \(e_i\)
\(\beta_i\) \(\sim\) N(Z\(\Delta\)[i,], \(V_{\beta}\))
Note: Z\(\Delta\) is the matrix \(Z * \Delta\) and [i,] refers to \(i\)th row of this product

\(vec(\Delta)\) given \(V_{\beta}\) \(\sim\) \(N(vec(Deltabar), V_{\beta}(x) A^{-1})\)
\(V_{\beta}\) \(\sim\) \(IW(nu,V)\)
\(Delta, Deltabar\) are \(nz x nvar\); \(A\) is \(nz x nz\); \(V_{\beta}\) is \(nvar x nvar\).

Note: if you don't have any \(Z\) variables, omit \(Z\) in the Data argument and a vector of ones will be inserted; the matrix \(\Delta\) will be \(1 x nvar\) and should be interpreted as the mean of all unit \(\beta\)s.

Argument Details

Data = list(regdata, Z) [Z optional]

regdata: list of lists with \(X\) and \(y\) matrices for each of nreg=length(regdata) regressions
regdata[[i]]$X: \(n_i x nvar\) design matrix for equation \(i\)
regdata[[i]]$y: \(n_i x 1\) vector of observations for equation \(i\)
Z: \(nreg x nz\) matrix of unit characteristics (def: vector of ones)

Prior = list(Deltabar, A, nu.e, ssq, nu, V) [optional]

Deltabar: \(nz x nvar\) matrix of prior means (def: 0)
A: \(nz x nz\) matrix for prior precision (def: 0.01I)
nu.e: d.f. parameter for regression error variance prior (def: 3)
ssq: scale parameter for regression error var prior (def: var(y_i))
nu: d.f. parameter for Vbeta prior (def: nvar+3)
V: Scale location matrix for Vbeta prior (def: nu*I)

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

R: number of MCMC draws
keep: MCMC thinning parm -- 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)

References

For further discussion, see Chapter 3, Bayesian Statistics and Marketing by Rossi, Allenby, and McCulloch.

See Also

rhierLinearMixture

Examples

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

nreg = 100
nobs = 100
nvar = 3
Vbeta = matrix(c(1, 0.5, 0, 0.5, 2, 0.7, 0, 0.7, 1), ncol=3)

Z = cbind(c(rep(1,nreg)), 3*runif(nreg))
Z[,2] = Z[,2] - mean(Z[,2])
nz = ncol(Z)
Delta = matrix(c(1,-1,2,0,1,0), ncol=2)
Delta = t(Delta) # first row of Delta is means of betas
Beta = matrix(rnorm(nreg*nvar),nrow=nreg)%*%chol(Vbeta) + Z%*%Delta

tau = 0.1
iota = c(rep(1,nobs))
regdata = NULL
for (reg in 1:nreg) { 
  X = cbind(iota, matrix(runif(nobs*(nvar-1)),ncol=(nvar-1)))
	y = X%*%Beta[reg,] + sqrt(tau)*rnorm(nobs)
	regdata[[reg]] = list(y=y, X=X) 
}

Data1 = list(regdata=regdata, Z=Z)
Mcmc1 = list(R=R, keep=1)

out = rhierLinearModel(Data=Data1, Mcmc=Mcmc1)

cat("Summary of Delta draws", fill=TRUE)
summary(out$Deltadraw, tvalues=as.vector(Delta))

cat("Summary of Vbeta draws", fill=TRUE)
summary(out$Vbetadraw, tvalues=as.vector(Vbeta[upper.tri(Vbeta,diag=TRUE)]))

## plotting examples
if(0){
  plot(out$betadraw)
  plot(out$Deltadraw)
}

Run the code above in your browser using DataLab