<jupyter_start><jupyter_text>Classical Approximate Bayesian Compuation<jupyter_code>import hypothesis
import matplotlib.pyplot as plt
import torch<jupyter_output><empty_output><jupyter_text>As in every Bayesian inference setting we have to define a prior $p(\theta)$, which is in our case $\mathcal{U}(-10, 10)$.<jupyter_code>from torch.distributions.uniform import Uniform

prior = Uniform(-10, 10)<jupyter_output><empty_output><jupyter_text>The simulation-based or likelihood-free setting also requires the definition of a forward model.<jupyter_code>from hypothesis.simulation import Simulator

class NormalSimulator(Simulator):
    
    def __init__(self):
        super(NormalSimulator, self).__init__()
        
    def forward(self, inputs):
        inputs = inputs.view(-1, 1)

        return torch.randn(inputs.size(0), 1) + inputs

simulator = NormalSimulator()<jupyter_output><empty_output><jupyter_text>Let's check how the joint $p(\theta,x)$ looks like.<jupyter_code>from hypothesis.visualization.util import make_square

n = torch.Size([100000])
inputs = prior.sample(n)
outputs = simulator(inputs)

plt.hist(outputs.numpy(), bins=50, density=True, color="black", histtype="step", lw=2)
plt.xlabel(r"$\theta$")
plt.ylabel("Posterior mass function")
make_square(plt.gca())
plt.show()<jupyter_output><empty_output><jupyter_text>Every Bayesian setting requires an observation $x$ (or a set of observations) to compute the posterior of:<jupyter_code>truth = prior.sample()
observation = simulator(truth)

print("Observation:", observation.item())<jupyter_output>Observation: 5.455133438110352<jupyter_text>And finally the posterior sampling procedure itself. The `acceptor` method takes two summary statistics (the observation and the sample being tested), and determines whether or not to accept the sample that has been used to generate `b`.<jupyter_code>%%time

from hypothesis.inference.abc import ApproximateBayesianComputation as ABC

def acceptor(a, b, epsilon=1):
    return (a - b).abs().item() < epsilon

def plot_posterior(posterior_samples):
    plt.title("Approximate posterior")
    plt.hist(posterior_samples.numpy(), bins=50, color="black", histtype="step", density=True, lw=2)
    plt.xlabel(r"$\theta$")
    plt.ylabel("Posterior mass function")
    plt.axvline(truth, c="C0", lw=2)
    make_square(plt.gca())
    plt.show()

# Use the identity function as a summary statistic.
def identity(x):
    return x
    
# Allocate the inference procedure.
abc = ABC(simulator=simulator,
    prior=prior,
    summary=identity,
    acceptor=acceptor)

# Sample the approximate posterior. We don't want gradients here (but you can).
with torch.no_grad():
    posterior_samples = abc.sample(observation, num_samples=10000)
    
# Plot the approximate posterior.
plot_posterior(posterior_samples)<jupyter_output><empty_output><jupyter_text>Let's check what the effect of a smaller `epsilon` value is in the acceptor method. This might take a bit longer.<jupyter_code>%%time

def acceptor(a, b, epsilon=0.05):
    return (a - b).abs().item() < epsilon

# Allocate the inference procedure with the new acceptor function.
abc = ABC(simulator=simulator,
    prior=prior,
    summary=identity, 
    acceptor=acceptor)

# Sample the approximate posterior. We don't want gradients here (but you can).
with torch.no_grad():
    posterior_samples = abc.sample(observation, num_samples=10000)
    
# Plot the approximate posterior.
plot_posterior(posterior_samples)<jupyter_output><empty_output><jupyter_text>The `ABC` procedure typically uses a single thread to sample from the posterior. However, the sampler can easily be parallelized.<jupyter_code>%%time

from hypothesis.inference.abc import ParallelApproximateBayesianComputation as ParallelABC

# Allocate the parallel ABC sampler.
parallel_abc_sampler = ParallelABC(abc, workers=8)

# Sample the approximate posterior. We don't want gradients here (but you can).
with torch.no_grad():
    posterior_samples = parallel_abc_sampler.sample(observation, num_samples=10000)
    
# Plot the approximate posterior.
plot_posterior(posterior_samples)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10350665054413542, 'fraction_numerical': 0.014993954050785973, 'mean_word_length': 5.692556634304207, 'pattern_counts': {'":': 0, '<': 34, '<?xml version=': 0, '>': 32, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '571444', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 1345, 'n_tokens_neox': 1250, 'n_words': 350, 'path': 'examples/approximate-bayesian-computation/classical_abc.ipynb', 'repo_name': 'boyali/hypothesis-sre', 'size': '4135'}