Understanding Bayesian Histogram

Question

I'm reading about Dirichlet Process Models in "Bayesian Data Analysis" by Gelman et. al.

To motivate the idea, they start with a section on Bayesian histograms. I am a little confused about their discussion. I'm sure it's obvious to the initiated, but alas, I'm uninitiated. I just want to verify that I understand correctly. I understand it might be hard to give an "answer" that's not just "yes" if I'm right. However, if enough people just say "yes", in the comments, I can answer my question by a simulated illustration, or something. Of course, if I'm wrong, then an answer would be correcting me :).

Summary of their description

We have iid data $y_i \sim f$, and we divide the data into intervals defined by $\xi_0 < \xi_1 < \cdots < \xi_k $. That is, $y_i \in [\xi_0,\xi_k]$.

The probability model for the density is

$$ f(y)=\sum_{h=1}^{k} 1_{\xi_{h-1}<y\le\xi_h}\frac{\pi_h}{(\xi_h-\xi_{h-1})}, y\in \mathbb{R}$$

Where $\pi=(\pi_1,\ldots,\pi_k)$ is the parameter of probabilities. We can work out the posterior of $\pi$ under a $Dirichlet(a_1 \ldots a_k)$ prior.

$$p(\pi \mid y) =^{\mathcal{D}} Dirichlet(a_1 + n_1, \ldots, a_k + n_k)$$

where $n_h = \sum_i 1_{\xi_{h-1} < y_i \le \xi_h} $

My Understanding

If we want to get the "Bayes histogram" density estimate, we

1. Specify the hyperparameters $a_1 \ldots a_k$
2. We calculate the values of $n_h$ from the data
3. We sample a bunch of $\pi$ vectors from the posterior distribution and calculate the posterior mean vecctor (or some other estimate)
4. We take $\hat{f}(y) = \frac{\pi_h}{(\xi_h-\xi_{h-1})}$ as our estimated density over the interval $[\xi_{h-1}, \xi_{h}]$, ie for all $y\in [\xi_{h-1}, \xi_{h}]$

Answer 1

Yes, I think your interpretation is correct. An important point is that the parameter vector $\pi$ in the model is assumed to follow a Dirichlet distribution. These parameters represents the probability that an observation ends up in each bin of the histogram, meaning that $\pi_h$ would be the probability of an observations $y_i$ being in the interval $[\xi_{h-1}, \xi_{h}]$. The hyperparameters $a_i$ can be seen as pseudo-counts, meaning how many observations we believe should be in the histogram bins before we seen any data.

The update is done as you described, and we get the posterior distribution for the parameter vector. The posterior distribution is also a Dirichlet distribution since the Dirichlet distribution is the conjugate prior for the likelihood function in the model. Then we can sample parameter vectors from the posterior, and do whatever analysis we want on that sample.

If you are reading the 3rd edition of "Bayesian Data Analysis" I suggest you check out the section 3.4 on the Multinomial model in chapter 3, which also uses a Dirichlet distribution as prior and works in a similar way.

Source: Bayesian Data Analysis 3rd edn A. Gelman, J. B. Carlin, H. S. Stern, D. B. Dunson, A. Vehtari and D. B. Rubin, 2013 Boca Raton, Chapman and Hall–CRC 676 pp., ISBN 1‐439‐84095‐4