# Bayesian Hierarchical Clustering: How to calculate probability of Data under $H_1$?

I am working on implementing the Bayesian hierarchical clustering algorithm found here from scratch as a way to practice and learn the algorithm. However, I have hit a snag in calculating the quantity in formula 1 of section 2:

\begin{align} p(D_k|H_1^k) &= \int p(D_k|\theta) p(\theta |\beta) d\theta\\ &= \int \left[ \prod_{x_i \in D_k} p(x_i |\theta) \right] p(\theta | \beta) d\theta \end{align}

In the subsequent paragraph, the authors state,

This is a natural model-based criterion for measuring how well the data fit into one cluster. If we choose models with conjugate priors (e.g. Normal-Inverse-Wishart priors for Normal continuous data or Dirichlet priors for Multinomial discrete data) this integral is tractable

I understand that this is the probability of cluster $$k$$ being generated from some distribution with parameter(s) $$\mathbf{\theta}$$ and that, on the first iteration of the algorithm, I should be calculating this quantity for every pair of data points (i.e. vectors) in the set.

However, I am a bit confused regarding the quote above. Does this mean that I should replace the contents of the integral with an appropriate posterior distribution and then perform the integration with respect to $$\mathbf{\theta}$$?

What would this mean in the Gaussian example mentioned by the authors:

For example, in the case of Gaussians, $$p(D_k|H_1^k)$$ is a function of the sample mean and covariance of the data in $$D_k$$

The paper seems to allude to this being a straight-forward calculation, but I have not seen conjugate priors used in this way. Any help would be appreciated!

I've been messing around with this a little bit myself, and what's being referenced here is the fact that the marginal likelihood (the integral expression you provide) has a quick-to-evaluate closed form that only makes use of the sample statistics of $$D_{k}$$ when the conjugate prior is employed, no numerical integration required.
One more tangentially-related thing. Practically speaking, you'll actually want to be dealing with the marginal log-likelihood, as the gamma function risks overflow. I heartily recommend you compute the (pseudo?)log-odds to make BHC merge decisions: $$\log(P(D_{k}|H_{1}^{k})P(H_{1}^{k})) - \log(P(D_{k}|H_{2}^{k})P(H_{2}^{k}))$$ Then, you can just plug the marginal log-likelihood right in.