For a binary state, its expected value is its probability of being one. Hence the Monte Carlo method estimate of the posterior methylation profile under the chosen target isThe averaging axis must be zero: average over sampled patterns while keeping the site index. A provisional choice for this particular trace is , retaining 2000 states:
def mean_methylation(meths, burnin):
if not 0 <= burnin < len(meths):
raise ValueError("burn-in must leave at least one state")
return np.asarray(meths)[burnin:].mean(axis=0)
burnin = 3000
profile_estimate = mean_methylation(meths, burnin)
print(lls[burnin:].mean())
plt.plot(profile_estimate, label="retained-sample mean")
plt.scatter(np.arange(len(obs)), .2*obs, s=10, label="0.2 × observation")
plt.xlabel("Site index")
plt.ylabel("Methylation probability")
plt.legend()
plt.show()The retained mean log score is , and the profile is shown below the trace in part (b). The fixed endpoints have estimate exactly one; adjacent sites are often smoothed into methylated or unmethylated stretches despite individual conflicting observations. The rescaled observations are only a visual comparison and are not another estimate of the posterior probability.
This is a finite-chain estimate for the conditional-product target defined in part (a). A burn-in chosen by inspection does not establish convergence, and the long excursions in this short Markov chain leave substantial uncertainty in individual-site estimates; several longer runs would be needed for precise inference. To sample the coherent joint alternative instead, pass
score=joint_methylation_score and inspect that new trace before choosing its own burn-in. Posterior methylation profile 2026-10-20
A posterior methylation profile gives the posterior probability of DNA methylation at each modeled site. It differs from one sampled binary pattern and from a maximum-probability joint pattern.