Computing the Bayesian Factor from a Markov chain Monte Carlo Simulation
  of the Posterior Distribution
  Computation of the marginal likelihood from a simulated posterior
distribution is central to Bayesian model selection but is computationally
difficult. I argue that the marginal likelihood can be reliably computed from a
posterior sample by careful attention to the numerics of the probability
integral. Posing the expression for the marginal likelihood as a Lebesgue
integral, we may convert the harmonic mean approximation from a sample
statistic to a quadrature rule. As a quadrature, the harmonic mean
approximation suffers from enormous truncation error as consequence . In
addition, I demonstrate that the integral expression for the harmonic-mean
approximation converges slowly at best for high-dimensional problems with
uninformative prior distributions. These observations lead to two
computationally-modest families of quadrature algorithms that use the full
generality sample posterior but without the instability. The first algorithm
automatically eliminates the part of the sample that contributes large
truncation error. The second algorithm uses the posterior sample to assign
probability to a partition of the sample space and performs the marginal
likelihood integral directly. This eliminates convergence issues. The first
algorithm is analogous to standard quadrature but can only be applied for
convergent problems. The second is a hybrid of cubature: it uses the posterior
to discover and tessellate the subset of that sample space was explored and
uses quantiles to compute a representive field value. Neither algorithm makes
strong assumptions about the shape of the posterior distribution and neither is
sensitive outliers. [abridged]
