Calculating 95% confidence intervals in quantile regression in R using rq function

Question

I would like to get 95% confidence intervals for the regression coefficients of a quantile regression. You can calculate quantile regressions using the rq function of the quantreg package in R (compared to an OLS model):

library(quantreg)
LM<-lm(mpg~disp, data = mtcars)
QR<-rq(mpg~disp, data = mtcars, tau=0.5)

I am able to get 95% confidence intervals for the linear model using the confint function:

confint(LM)

When I use quantile regression I understand that the following code produces bootstrapped standard errors:

summary.rq(QR,se="boot")

But actually I would like something like 95% confidence intervals. That is, something to interprete like: "with a probability of 95%, the interval [...] includes the true coefficient". When I calculate standard errors using summary.lm() I would just multiply SE*1.96 and get similar results as from confint(). But this is not possible using bootstrapped standard errors. So my question is how get 95% confidence intervals for quantile regression coefficients?

How wrong is it to compute a 95% CI this way: estimate +- 1.96 * SE ? The software reports the SEs and pvalues, so presumably the SEs can be used to form CIs as well. — frelk, Dec 21 '21 at 17:06

score 8 · Accepted Answer · answered Jun 29 '16 at 18:06

You can use the boot.rq function directly to bootstrap the coefficients:

x<-1:50
y<-c(x[1:48]+rnorm(48,0,5),rnorm(2,150,5))

QR <- rq(y~x, tau=0.5)
summary(QR, se='boot')

LM<-lm(y~x)

QR.b <- boot.rq(cbind(1,x),y,tau=0.5, R=10000)

t(apply(QR.b$B, 2, quantile, c(0.025,0.975)))
confint(LM)


plot(x,y)
abline(coefficients(LM),col="green")
abline(coefficients(QR),col="blue")

for(i in seq_len(nrow(QR.b$B))) {
  abline(QR.b$B[i,1], QR.b$B[i,2], col='#0000ff01')
}

You may want to use the boot package to compute intervals other than the percentile interval.

Nice answer. In the t(apply...) `Error in QR.b$B : $ operator is invalid for atomic vectors`. I think no column selection for the QR.b: `t(apply(QR.b, 2, quantile, c(0.025,0.975)))`? — Minnow, Apr 24 '17 at 15:30

score 1 · Answer 2 · answered Nov 16 '17 at 00:44

You could also simply retrieve the vcov from the object, setting covariance=TRUE. This amounts to using boostrapped standard errors in your CI:

vcov.rq <- function(x, se = "iid") {
 vc <- summary.rq(x, se=se, cov=TRUE)$cov
 dimnames(vc) <- list(names(coef(x)), names(coef(x)))
 vc
}

confint(QR)

But yes, the better way to do this is to use a studentized bootstrap.

Calculating 95% confidence intervals in quantile regression in R using rq function

2 Answers2

Linked