# \donttest{
require(spatialprobit)
# number of observations
n <- 100
# true parameters
beta <- c(0, 1, -1)
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!
W <- kNearestNeighbors(x = rnorm(n), y = rnorm(n), k = m)
# innovations
eps <- rnorm(n=n, mean=0, sd=1)
# generate data from model
S <- I_n - rho * W
z <- solve(qr(S), X %*% beta + eps)
y <- as.vector(z >= 0) # 0 or 1, FALSE or TRUE
# estimate SAR probit model
set.seed(12345)
sarprobit.fit1 <- sar_probit_mcmc(y, X, W, ndraw=500, burn.in=100,
thinning=1, prior=NULL, computeMarginalEffects=TRUE)
summary(sarprobit.fit1)
# print impacts
impacts(sarprobit.fit1)
################################################################################
#
# Example from LeSage/Pace (2009), section 10.3.1, p. 302-304
#
################################################################################
# Value of "a" is not stated in book!
# Assuming a=-1 which gives approx. 50% censoring
library(spatialprobit)
a <- -1 # control degree of censored observation
n <- 1000
rho <- 0.7
beta <- c(0, 2)
sige <- 0.5
I_n <- sparseMatrix(i=1:n, j=1:n, x=1)
x <- runif(n, a, 1)
X <- cbind(1, x)
eps <- rnorm(n, sd=sqrt(sige))
param <- c(beta, sige, rho)
# random locational coordinates and 6 nearest neighbors
lat <- rnorm(n)
long <- rnorm(n)
W <- kNearestNeighbors(lat, long, k=6)
y <- as.double(solve(I_n - rho * W) %*% (X %*% beta + eps))
table(y > 0)
# set negative values to zero to reflect sample truncation
ind <- which(y <=0)
y[ind] <- 0
# Fit SAR Tobit (with approx. 50% censored observations)
fit_sartobit <- sartobit(y ~ x, W, ndraw=1000, burn.in=200,
computeMarginalEffects=TRUE, showProgress=TRUE)
# print impacts
impacts(fit_sartobit)
# }
Run the code above in your browser using DataLab