CODEWITHSUNDEEP
rggplot2regression
Petr
Petr

Reputation: 1817

Plotting quantile regression sensitivity in ggplot

I would like to ask for advice how to plot quantile standart errors like with basic function in quantreg package in R

library(quantreg)
library(ggplot2)
library(dplyr)

QR.2 <- rq(hp ~  disp + mpg + I(mpg^2) + qsec + am, data = mtcars, tau = 1:9/10)
plot(summary(QR.2, se="boot"), ols=T)

Plots above shows quantile standart errors and confidence intervals. Is there a way how to reproduce same plot in ggplot?

Trying code below is pretty good start however geom_ribbon(aes(ymin=conf.low,ymax=conf.high),alpha=0.25, fill="#27408b") does return "confidence intervals" however clearly those are not the same as produced in plot above.

Is there a way how to get conf. intervals as in plot above?

rq(data=mtcars,
   tau= 1:9/10,
   formula = hp ~  disp + mpg + I(mpg^2) + qsec + am) %>% 
  broom::tidy() %>% 
  filter(!grepl("factor", term)) %>%   
  filter(!grepl("Intercept", term)) %>%   
  ggplot(aes(x=tau,y=estimate))+
  geom_point(color="#27408b", size = 3)+ 
  geom_line(color="#27408b", size = 1)+ 
  geom_smooth(method=  "lm", colour = "red", se = T)+  
  facet_wrap(~term, scales="free", ncol=2) + 
  geom_ribbon(aes(ymin=conf.low,ymax=conf.high),alpha=0.25, fill="#27408b")

Upvotes: 1

Views: 1765

Answers (3)

Gonzalo Rizzo
Gonzalo Rizzo

Reputation: 151

The accepted answer does not address the asked question since it does not plot the OLS point estimate (red) with ci (dashed red), as it is done by quantreg package (as Perl highlighted). The correct answer to this question would be the following:

# OLS
lm <- lm(data=mtcars, 
         formula =  hp ~  disp + mpg + I(mpg^2) + qsec + am)

ols <- as.data.frame(coef(lm))
ols.ci <- as.data.frame(confint(lm))
ols2 <- cbind(ols, ols.ci)
ols2 <- tibble::rownames_to_column(ols2, var="term")

# Quantile
rq(data=mtcars, 
   tau= 1:9/10,
   formula = hp ~  disp + mpg + I(mpg^2) + qsec + am) %>%
  broom::tidy(se.type = "boot") %>%
  filter(!grepl("factor", term)) %>%
  ggplot(aes(x=tau,y=estimate))+
  
  # quantilie results
  geom_point(color="#27408b", size = 3)+ 
  geom_line(color="#27408b", size = 1)+ 
  geom_ribbon(aes(ymin=conf.low,ymax=conf.high),alpha=0.25, fill="#27408b")+

  # OLS results
  geom_hline(data = ols2, aes(yintercept= `coef(lm)`), lty=1, color="red", size=1)+
  geom_hline(data = ols2, aes(yintercept= `2.5 %`), lty=2, color="red", size=1)+
  geom_hline(data = ols2, aes(yintercept= `97.5 %`), lty=2, color="red", size=1)+
  facet_wrap(~term, scales="free", ncol=2)

Upvotes: 1

Harry Haupt
Harry Haupt

Reputation: 21

The plot of the nine quantile regression models estimated by Petr displays, for each constant marginal effect, the qr point estimate (for each of the selected nine deciles) and a corresponding confidence interval (pointwise, with default alpha). Petr has selected to also display the OLS point estimate (red), again with ci (dashed red), a regular practice often used to compare the conditional median to the conditional mean. The latter is displayed correctly - as a horizontal line - in the plot generated by quantreg. The red line displayed using ggplot2 in the answer is obviously something different.

Upvotes: 1

Jon Spring
Jon Spring

Reputation: 66540

To get the equivalent, I think you need to use broom::tidy(se.type = "boot") %>%, otherwise the standard errors are calculated using a different method.

Base R output: enter image description here

ggplot2 equiv: enter image description here

  rq(data=mtcars, 
     tau= 1:9/10,
     formula = hp ~  disp + mpg + I(mpg^2) + qsec + am) %>%
  broom::tidy(se.type = "boot") %>%
  filter(!grepl("factor", term)) %>%
  ggplot(aes(x=tau,y=estimate))+
  geom_point(color="#27408b", size = 3)+ 
  geom_line(color="#27408b", size = 1)+ 
  geom_smooth(method=  "lm", colour = "red", se = T)+  
  facet_wrap(~term, scales="free", ncol=2) + 
  geom_ribbon(aes(ymin=conf.low,ymax=conf.high),alpha=0.25, fill="#27408b")

Upvotes: 3

Related Questions