HeckmanStan {HeckmanStan} | R Documentation |
Fit the Heckman Selection Stan model using the Normal, Student-t or Contaminated Normal distributions.
Description
'HeckmanStan()' fits the Heckman selection model using a Bayesian approach to address sample selection bias.
Usage
HeckmanStan(
y,
x,
w,
cc,
family = "CN",
init = "random",
thin = 5,
chains = 1,
iter = 10,
warmup = 5
)
Arguments
y |
A response vector. |
x |
A covariate matrix for the response y. |
w |
A covariate matrix for the missing indicator cc. |
cc |
A missing indicator vector (1=observed, 0=missing) . |
family |
The distribution family to be used (Normal, T, or CN). |
init |
Parameters specifies the initial values for model parameters. |
thin |
An Interval at which samples are retained from the MCMC process to reduce autocorrelation. |
chains |
The number of chains to run during the MCMC sampling. Running multiple chains is useful for checking convergence. |
iter |
The total number of iterations for the MCMC sampling, determining how many samples will be drawn. |
warmup |
The number of initial iterations that will be discarded as the algorithm stabilizes before collecting samples. |
Value
An object of class HeckmanStan
, which is a list containing two elements:
-
list[[1]]
: Includes inference results from the Stan model, along with EAIC and EBIC. -
list[[2]]
: Includes the HPC confidence intervals, along with LOOIC, WAIC, and CPO.
Examples
################################################################################
# Simulation
################################################################################
library(mvtnorm)
n<- 100
w<- cbind(1,rnorm(n),rnorm(n))
x<- cbind(w[,1:2])
family="CN"
sigma2<- 1
rho<-0.7
beta<- c(1,0.5)
gamma<- c(1,0.3,-.5)
nu=c(0.1,0.1)
data<-geraHeckman(x,w,beta,gamma,sigma2,rho,nu,family=family)
y<-data$y
cc<-data$cc
# Fit Heckman Normal Stan model
fit.n_stan <- HeckmanStan(y, x, w, cc, family="Normal"
, thin = 5, chains = 1, iter = 10000, warmup = 1000)
qoi=c("beta","gamma","sigma_e","sigma2", "rho","EAIC","EBIC")
print(fit.n_stan[[1]],par=qoi)
print(fit.n_stan[[2]])
require(rstan)
plot(fit.n_stan[[1]], pars=qoi)
plot(fit.n_stan[[1]], plotfun="hist", pars=qoi)
plot(fit.n_stan[[1]], plotfun="trace", pars=qoi)
plot(fit.n_stan[[1]], plotfun = "rhat")