Implement sampling from posterior + prediction intervals

If I didn't overlook something, MCMC currently only gives you uncertainties in the parameters? Functionality to sample from posterior + plot prediction intervals would be good. emcee seems to implement posterior sampling like this:

flat_samples = sampler.get_chain(discard=100, thin=15, flat=True)
inds = np.random.randint(len(flat_samples), size=100)
for ind in inds:
    sample = flat_samples[ind]
    # do something with sample here

https://emcee.readthedocs.io/en/stable/tutorials/line/