bayesm (version 3.1-6)

rmultireg: Draw from the Posterior of a Multivariate Regression

Description

rmultireg draws from the posterior of a Multivariate Regression model with a natural conjugate prior.

Usage

rmultireg(Y, X, Bbar, A, nu, V)

Value

A list of the components of a draw from the posterior

B

draw of regression coefficient matrix

Sigma

draw of Sigma

Arguments

Y

\(n x m\) matrix of observations on m dep vars

X

\(n x k\) matrix of observations on indep vars (supply intercept)

Bbar

\(k x m\) matrix of prior mean of regression coefficients

A

\(k x k\) Prior precision matrix

nu

d.f. parameter for Sigma

V

\(m x m\) pdf location parameter for prior on Sigma

Warning

This routine is a utility routine that does not check the input arguments for proper dimensions and type.

Author

Peter Rossi, Anderson School, UCLA, .

Details

Model:
\(Y = XB + U\) with \(cov(u_i) = \Sigma\)
\(B\) is \(k x m\) matrix of coefficients; \(\Sigma\) is \(m x m\) covariance matrix.

Priors:
\(\beta\) | \(\Sigma\) \(\sim\) \(N(betabar, \Sigma(x) A^{-1})\)
\(betabar = vec(Bbar)\); \(\beta = vec(B)\)
\(\Sigma\) \(\sim\) IW(nu, V)

References

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

Examples

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

n =200
m = 2
X = cbind(rep(1,n),runif(n))
k = ncol(X)
B = matrix(c(1,2,-1,3), ncol=m)
Sigma = matrix(c(1, 0.5, 0.5, 1), ncol=m)
RSigma = chol(Sigma)
Y = X%*%B + matrix(rnorm(m*n),ncol=m)%*%RSigma

betabar = rep(0,k*m)
Bbar = matrix(betabar, ncol=m)
A = diag(rep(0.01,k))
nu = 3
V = nu*diag(m)

betadraw = matrix(double(R*k*m), ncol=k*m)
Sigmadraw = matrix(double(R*m*m), ncol=m*m)

for (rep in 1:R) {
  out = rmultireg(Y, X, Bbar, A, nu, V)
  betadraw[rep,] = out$B
  Sigmadraw[rep,] = out$Sigma
  }

cat(" Betadraws ", fill=TRUE)
mat = apply(betadraw, 2, quantile, probs=c(0.01, 0.05, 0.5, 0.95, 0.99))
mat = rbind(as.vector(B),mat)
rownames(mat)[1] = "beta"
print(mat)

cat(" Sigma draws", fill=TRUE)
mat = apply(Sigmadraw, 2 ,quantile, probs=c(0.01, 0.05, 0.5, 0.95, 0.99))
mat = rbind(as.vector(Sigma),mat); rownames(mat)[1]="Sigma"
print(mat)