semprobit function

Bayesian estimation of the SEM probit model

Bayesian estimation of the probit model with spatial errors (SEM probit model).

semprobit(formula, W, data, subset, ...)

sem_probit_mcmc(y, X, W, ndraw = 1000, burn.in = 100, thinning = 1, 
  prior=list(a1=1, a2=1, c=rep(0, ncol(X)), T=diag(ncol(X))*1e12,
  nu = 0, d0 = 0, lflag = 0), 
  start = list(rho = 0.75, beta = rep(0, ncol(X)), sige = 1),
  m=10, showProgress=FALSE, univariateConditionals = TRUE)

Arguments

• y: dependent variables. vector of zeros and ones

• X: design matrix

• W: spatial weight matrix

• ndraw: number of MCMC iterations

• burn.in: number of MCMC burn-in to be discarded

• thinning: MCMC thinning factor, defaults to 1.

• prior: A list of prior settings for $or rho ~ Beta(a1,a2)$

  and $beta ~ N(c,T)$. Defaults to diffuse prior for beta.

• start: list of start values

• m: Number of burn-in samples in innermost Gibbs sampler. Defaults to 10.

• showProgress: Flag if progress bar should be shown. Defaults to FALSE.

• univariateConditionals: Switch whether to draw from univariate or multivariate truncated normals. See notes. Defaults to TRUE.

• formula: an object of class "formula" (or one that can be coerced to that class): a symbolic description of the model to be fitted.

• data: an optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the variables in the model. If not found in data, the variables are taken from environment(formula), typically the environment from which semprobit is called.

• subset: an optional vector specifying a subset of observations to be used in the fitting process.

• ...: additional arguments to be passed

Details

Bayesian estimates of the probit model with spatial errors (SEM probit model)

which leads to the data-generating process

where y is a binary 0,1 $(n x 1)$ vector of observations for $z < 0$ and $z >= 0$. $beta$ is a $(k x 1)$

vector of parameters associated with the $(n x k)$ data matrix X.

The prior distributions are $beta ~ N(c,T)$, $\sigma^2_{\epsilon} \sim IG(a1, a2)$, and $rho ~ Uni(rmin,rmax)$

or $or rho ~ Beta(a1,a2)$.

Returns

Returns a structure of class semprobit:

• beta: posterior mean of bhat based on draws

• rho: posterior mean of rho based on draws

• bdraw: beta draws (ndraw-nomit x nvar)

• pdraw: rho draws (ndraw-nomit x 1)

• sdraw: sige draws (ndraw-nomit x 1)

• total: a matrix (ndraw,nvars-1) total x-impacts

• direct: a matrix (ndraw,nvars-1) direct x-impacts

• indirect: a matrix (ndraw,nvars-1) indirect x-impacts

• rdraw: r draws (ndraw-nomit x 1) (if m,k input)

• nobs: # of observations

• nvar: # of variables in x-matrix

• ndraw: # of draws

• nomit: # of initial draws omitted

• nsteps: # of samples used by Gibbs sampler for TMVN

• y: y-vector from input (nobs x 1)

• zip: # of zero y-values

• a1: a1 parameter for beta prior on rho from input, or default value

• a2: a2 parameter for beta prior on rho from input, or default value

• time: total time taken

• rmax: 1/max eigenvalue of W (or rmax if input)

• rmin: 1/min eigenvalue of W (or rmin if input)

• tflag: 'plevel' (default) for printing p-levels; 'tstat' for printing bogus t-statistics

• lflag: lflag from input

• cflag: 1 for intercept term, 0 for no intercept term

• lndet: a matrix containing log-determinant information (for use in later function calls to save time)

References

LeSage, J. and Pace, R. K. (2009), Introduction to Spatial Econometrics, CRC Press, chapter 10

Author(s)

adapted to and optimized for R by Stefan Wilhelm wilhelm@financial.com based on code from James P. LeSage

See Also

sar_lndet for computing log-determinants

Examples

library(Matrix)
# number of observations
n <- 200

# true parameters
beta <- c(0, 1, -1)
sige <- 2
rho <- 0.75

# design matrix with two standard normal variates as "covariates"
X <- cbind(intercept=1, x=rnorm(n), y=rnorm(n))

# sparse identity matrix
I_n <- sparseMatrix(i=1:n, j=1:n, x=1)

# number of nearest neighbors in spatial weight matrix W
m <- 6

# spatial weight matrix with m=6 nearest neighbors
# W must not have non-zeros in the main diagonal!
i <- rep(1:n, each=m)
j <- rep(NA, n * m)
for (k in 1:n) {
  j[(((k-1)*m)+1):(k*m)] <- sample(x=(1:n)[-k], size=m, replace=FALSE)
}
W <- sparseMatrix(i, j, x=1/m, dims=c(n, n))

# innovations
eps <- sqrt(sige)*rnorm(n=n, mean=0, sd=1)

# generate data from model 
S <- I_n - rho * W
z <- X %*% beta + solve(qr(S), eps)
y <- as.double(z >= 0)  # 0 or 1, FALSE or TRUE

# estimate SEM probit model
semprobit.fit1 <- semprobit(y ~ X - 1, W, ndraw=500, burn.in=100, 
  thinning=1, prior=NULL)
summary(semprobit.fit1)