rsurGibbs {bayesm} | R Documentation

## Gibbs Sampler for Seemingly Unrelated Regressions (SUR)

### Description

`rsurGibbs` implements a Gibbs Sampler to draw from the posterior of the Seemingly Unrelated Regression (SUR) Model of Zellner

### Usage

    rsurGibbs(Data, Prior, Mcmc)

### Arguments

`Data` | list(regdata)
`Prior` | list(betabar,A, nu, V)
`Mcmc` | list(R,keep)

### Details

Model: y_i = X_ibeta_i + e_i. i=1,...,m. m regressions.
(e(1,k), ..., e(m,k)) ~ N(0,Sigma). k=1, ..., nobs.

We can also write as the stacked model:
y = Xbeta + e where y is a nobs*m long vector and k=length(beta)=sum(length(betai)).

Note: we must have the same number of observations in each equation but we can have different numbers of X variables

Priors: beta ~ N(betabar,A^{-1}). Sigma ~ IW(nu,V).

List arguments contain

`regdata`
  list of lists, regdata[[i]]=list(y=yi,X=Xi)
`betabar`
  k x 1 prior mean (def: 0)
`A`
  k x k prior precision matrix (def: .01I)
`nu`
  d.f. parm for Inverted Wishart prior (def: m+3)
`V`
  scale parm for Inverted Wishart prior (def: nu*I)
`R`
  number of MCMC draws
`keep`
  thinning parameter - keep every keepth draw

### Value

list of MCMC draws

`betadraw` | R x k array of betadraws
`Sigmadraw` | R x (m*m) array of Sigma draws

### Author(s)

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 and McCulloch, Chapter 3.
http://gsbwww.uchicago.edu/fac/peter.rossi/research/bsm.html

### See Also

`rmultireg`

### Examples

    if(nchar(Sys.getenv("LONG_TEST")) != 0) {R=1000} else {R=10}
    ##
    ## simulate data from SUR
    set.seed(66)
    beta1=c(1,2)
    beta2=c(1,-1,-2)
    nobs=100
    nreg=2
    iota=c(rep(1,nobs))
    X1=cbind(iota,runif(nobs))
    X2=cbind(iota,runif(nobs),runif(nobs))
    Sigma=matrix(c(.5,.2,.2,.5),ncol=2)
    U=chol(Sigma)
    E=matrix(rnorm(2*nobs),ncol=2)
    y1=X1%*%beta1+E[,1]
    y2=X2%*%beta2+E[,2]
    ##
    ## run Gibbs Sampler
    regdata=NULL
    regdata[[1]]=list(y=y1,X=X1)
    regdata[[2]]=list(y=y2,X=X2)
    Mcmc=list(R=R)
    out=rsurGibbs(Data=list(regdata=regdata),Mcmc=Mcmc)
    ##
    ## summarize GS output
    if(R < 100) {begin=1} else {begin=100}
    end=Mcmc$R
    cat(" betadraws ",fill=TRUE)
    mat=apply(out$betadraw[begin:end,],2,quantile,probs=c(.01,.05,.5,.95,.99))
    mat=rbind(c(beta1,beta2),mat); rownames(mat)[1]="beta"; print(mat)
    cat(" Sigmadraws ",fill=TRUE)
    mat=apply(out$Sigmadraw[begin:end,],2,quantile,probs=c(.01,.05,.5,.95,.99))
    mat=rbind(as.vector(Sigma),mat); rownames(mat)[1]="sigma"; print(mat)

[Package bayesm version 2.0-8 Index]