gistfile1.txt

extract_sims <- function(modelInterval){
  # Grabs simulation draws from a predictInterval object.
  attr(modelInterval, "sim.results")
}

bind_sims <- function(sims, newdata){
  # Where you have simulated predictions based the model or new dataset, it can
  # be helpful to combine these with the dataset. 
  # There may also be multiple draws from the simulation, stored in separate
  # columns. These are melted into long form.
  b <- bind_cols(as.data.frame(sims), newdata) 
  melt(b, id.vars=names(newdata)) %>% rename(fitted = value, sim=variable) %>% 
    mutate(sim=as.numeric(factor(sim)))
}

preds <- predictInterval(model, 
                  newdata=newdata, 
                  n.sims=NSIMS, 
                  returnSims=T, 
                  include.resid.var=F)




preds.long <- bind_sims(extract_sims(rf2.pred.withresid), newdata)

# then you can plot or whatever...

ggplot(preds.long, aes(x,y)) ...