QRc {bayesQR}
Description
QRc implements a Bayesian method for estimating quantile regression parameters. This method uses the asymmetric Laplace distribution for this goal (see references). To improve the speed of the routine, the Markov Chain Monte Carlo (MCMC) part of the algorithm is programmed in Fortran and is called from within the R function QRc.
Usage
QRc(Data, Prior, Mcmc)
Arguments
- Data
- list(y, X, p)
- Prior
- list(beta0, V0, shape0, scale0)
- Mcmc
- list(R, keep)
Details
| Model: | y = Xbeta + e |
| e ~ ALD(location=0, scale=sigma, quantile=p) | |
| Priors: | beta ~ N(beta0, V0) |
| sigma ~ invGamma(shape0,scale0) |
List arguments contain:
- X : n x nvar(X) matrix of predictors (first column should be a vector of ones if intercept is desired)
- y : n x 1 vector of observations (dependent variable)
- p : quantile of interest (should be between 0 and 1)
- beta0 : nvar(X) x 1 vector of prior means (default: 0)
- V0 : nvar(X) x nvar(X) prior covariance matrix (default: 100*diag(ncol(X)))
- shape0 : shape parameter for inverse Gamma prior for sigma (default: 0.01)
- scale0 : scale parameter for inverse Gamma prior for sigma (default: 0.01)
- R : number of MCMC draws
- keep : thinning parameter, i.e. keep every keepth draw (default: 1)
Values
A list containing:
- method
- a string containing 'QRc'
- p
- the quantile that was estimated
- betadraw
- R/keep x nvar(X) matrix of beta draws
- sigmadraw
- R/keep vector of sigma draws
References
The algorithm is an improved (Gibbs instead of Metropolis-Hastings) implementation of:
Yu, K. and Moyeed, R. (2001). Bayesian quantile regression, Statistics and Probability Letters, 54, 437-447.
Also see:
Kozumi, H. and Kobayashi, G. (2011). Gibbs sampling methods for Bayesian quantile regression, Journal of Statistical Computation and Simulation, 81(11), 1565-1578.
Examples
# Simulate data from heteroskedastic regression n <- 200 X <- runif(n=n,min=0,max=10) X <- cbind(1,X) y <- 1 + 2*X[,2] + rnorm(n=n, mean=0, sd=.6*X[,2]) # Initiate plot ## Plot datapoints plot(X[,2], y, main="", cex=.6, xlab="X") # Initialize the inputs for QRc Data = list(y=y, X=X, p=.05) Prior = list(beta0=rep(0,ncol(X)),V0=100*diag(ncol(X)), shape0=0.01, scale0=0.01) Mcmc = list(R=5000, keep=1) # Write loop to analyze 5 quantiles for (i in 1:5) { if (i==1) p = .05 if (i==2) p = .25 if (i==3) p = .50 if (i==4) p = .75 if (i==5) p = .95 Data$p = p out = QRc(Data=Data, Prior=Prior, Mcmc=Mcmc) ## Add quantile regression lines to the plot (exclude first 500 burn-in draws) abline(a=mean(out$betadraw[500:Mcmc$R,1]),b=mean(out$betadraw[500:Mcmc$R,2]),lty=i,col=i) } # Estimate and plot OLS model outOLS = lm(y~X[,2]) abline(outOLS,lty=1,lwd=2,col=6) # Add legend to plot legend(x=0,y=max(y),legend=c(.05,.25,.50,.75,.95,"OLS"),lty=c(1,2,3,4,5,1),lwd=c(1,1,1,1,1,2),col=c(1:6),title="Quantile")
Documentation reproduced from package bayesQR, version 2.0. License: GPL (>= 2)
