Calcolo degli intervalli di confidenza del 95% nella regressione quantile in R utilizzando la funzione rq

Question

Vorrei ottenere intervalli di confidenza del 95% per i coefficienti di regressione di una regressione quantile. È possibile calcolare regressioni quantili utilizzando la funzione rq del pacchetto quantreg in R (rispetto ad un modello OLS):

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

io sono in grado di ottenere gli intervalli di confidenza al 95% per il modello lineare usando la funzione confint:

confint(LM) 

Quando uso la regressione quantile capisco che il seguente codice produce bootstrap degli errori standard:

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

Ma in realtà vorrei qualcosa come gli intervalli di confidenza al 95%. Cioè, qualcosa da interpretare come: "con una probabilità del 95%, l'intervallo [...] include il vero coefficiente". Quando calcolo gli errori standard usando summary.lm() vorrei solo moltiplicare SE * 1.96 e ottenere risultati simili a partire da confint(). Ma questo non è possibile usando errori standard bootstrap. Quindi la mia domanda è: come ottenere intervalli di confidenza al 95% per i coefficienti di regressione quantile?

Answer 1

È possibile utilizzare la funzione boot.rq direttamente bootstrap coefficienti:

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') 
} 

È possibile utilizzare il pacchetto di avvio per calcolare intervalli diversi dell'intervallo percentile.

Answer 2

Si potrebbe anche semplicemente recuperare il vcov dall'oggetto, impostando covariance=TRUE. Ciò equivale a utilizzare errori standard boostrapped nella vostra 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) 

Ma sì, il modo migliore per farlo è quello di utilizzare un bootstrap studentizzato.