Confidence interval bands in ggplot2 when using stat_quantile?

Question

I would like to add the median spline and corresponding confidence interval bands to a ggplot2 scatter plot. I am using the 'quantreg'-package, more specifically the rqss function (Additive Quantile Regression Smoothing).

In ggplot2 I am able to add the median spline, but not the confidence interval bands:

fig = ggplot(dd, aes(y = MeanEst, x = N, colour = factor(polarization)))
fig + stat_quantile(quantiles=0.5, formula = y ~ qss(x), method = "rqss") +
  geom_point()

The quantreg-package comes with its own plot function; plot.rqss. Where I am able to add the confidence bands (bands=TRUE):

plot(1, type="n", xlab="", ylab="", xlim=c(2, 12), ylim=c(-3, 0)) # empty plot
plotfigs = function(df) {
  rqss_model = rqss(df$MeanEst ~ qss(df$N))
  plot(rqss_model, bands=TRUE, add=TRUE, rug=FALSE, jit=FALSE)
  return(NULL)
}
figures = lapply(split(dd, as.factor(dd$polarization)), plotfigs)

However plot function that comes with the quantreg-package is not very flexible/well suited for my needs. Is it possible to get the confidence bands in a ggplot2 plot? Perhaps by mimicking the method used in the quantreg-package, or simply copying them from the plot?

Data: pastebin.

Answer 1

You almost have it. When you call

 plot(rqss_model, bands=TRUE, add=TRUE, rug=FALSE, jit=FALSE)

The function very helpfully returns the plotted data. All we do is grab the data frame. First a minor tweak to your function, return the data in a sensible way

plotfigs = function(df) {
  rqss_model = rqss(df$MeanEst ~ qss(df$N))
  band = plot(rqss_model, bands=TRUE, add=TRUE, rug=FALSE, jit=FALSE)
  data.frame(x=band[[1]]$x, low=band[[1]]$blo, high=band[[1]]$bhi, 
             pol=unique(df$polarization))
}

Next call the function and condense

figures = lapply(split(dd, as.factor(dd$polarization)), plotfigs)
bands = Reduce("rbind", figures)

Then use geom_ribbon to plot

## We inherit y and color, so have to set them to NULL
fig + geom_ribbon(data=bands, 
                  aes(x=x, ymin=low, ymax=high, 
                  y=NULL, color=NULL, group=factor(pol)), 
                  alpha=0.3)