Bayesian rejection sampling for a segregating-site count

With proper prior $\pi (\theta )$ independent of the neutral Kingman's coalescent times, propose $(\theta ,T)$ from their joint prior and accept with probability $e^{-\theta L/2}(\theta L/2{)}^k/k!$. The accepted density is proportional to the prior times the observed-count likelihood, hence is the posterior given $S=k$. The acceptance rate is the prior predictive probability $P(S=k)$. For $k>0$, divide the acceptance function by its maximum $M_k=e^{-k}{k}^k/k!$, achieved at $\theta L/2=k$, to improve the rate to $P(S=k)/M_k$ without altering the accepted distribution. For $k=0$, $M_0=1$.

For a proper prior distribution $\pi (\theta )$ independent of coalescent times $T$, let $q(T)={\prod }_{j=2}^n{\lambda }_j{e}^{-{\lambda }_j{T}_j}$ with ${\lambda }_j=\left(\genfrac{}{}{0ex}{}{j}{2}\right)$ and $L=\sum j{T}_j$. Observing $S=k$ under an infinite sites mutation model gives posterior density $\pi (\theta )q(T)e^{-\theta L/2}(\theta L/2{)}^k/(k!m_k)$, where $m_k$ is the prior predictive probability of $S=k$. Multiplication by the conditional Poisson mass and normalization is an application of Bayes' theorem.

// Target: statistical-inference.bigb